Hermitian and non-Hermitian covariance estimators for multivariate Gaussian and 
non-Gaussian assets from random matrix theory 



Andrzej Jarosj^ 



(N 

o 

o3 



O 

(N 

H 

d 

> 

00 

OV 

(N 

d 



X 



The Henryk Niewodniczanski Institute of Nuclear Physics, 
Polish Academy of Sciences, Radzikowskiego 152, 31-342 Krakow, Poland 

The random matrix theory method of planar Gaussian diagrammatic expansion is applied to 
find the mean spectral density of the Hermitian equal-time and non-Hermitian time-lagged cross- 
covariance estimators, firstly in the form of master equations for the most general multivariate 
Gaussian system, secondly for seven particular toy models of the true covariance function. For 
the simplest one of these models, the existing result is shown to be incorrect and the right one 
is presented, moreover its generalizations are accomplished to the exponentially-weighted moving 
average estimator as well as two non-Gaussian distributions, Student t and free Levy. The paper 
revolves around applications to financial complex systems, and the results constitute a sensitive 
probe of the true correlations present there. 
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I. INTRODUCTION 

An important problem in various fields of science fea- 
turing systems of a large number of time-dependent ob- 
jects is to unravel correlations between these objects, ei- 
ther at the same or at distinct time moments fSec. II A Tj) . 
This knowledge may be gained from the past realized 
values of the correlations, although it is obscured by a 
measurement noise caused by an insufficient amount of 
information in any such time series (Sec. II A 2|) . This 
paper discusses how to account for this noise in ap- 
plications to models of complex systems which exhibit 
properties such as sectors of increased cross-correlations 
(Sec. II B 1[) or vector autoregression type of temporal dy- 
namics (Sec. IIB"2|) . 

Often in the below discussion the reference point will 
be financial markets — systems which are very complex, 
relevant for everyday life, accessible to empirical study 
and with the underlying mechanisms full of question 
marks — however, the results are beneficial for all kinds 
of complex systems, such as in brain science [l[ or mete- 
orology 0. 



A. Gaussian covariance functions and their 
historical estimation 

1. Gaussian covariance functions 

Random matrix theory approach to complex systems. 
Consider a system of N real variables (labeled by i = 
1,2, ... , N), and each one of them is observed through- 
out consecutive time moments a = 1, 2, . . . , T (separated 



by some time St); let the value of object i at time a be 
denoted by Ri a , all of which constitute an N x T ma- 
trix R. For such a large complex system, a standard 
approach is to replace complexity by randomness, i.e., 
assume that R is a random matrix, endowed with some 
joint probability density function (JPDF) V(R), and in- 
vestigate its statistical properties (here: the mean density 
of its eigenvalues), which then should help understand 
some of the generic features of the underlying complex 
system. Since there are two discrete dimensions here — 
spatial and temporal — random matrix theory (RMT) is 
a natural framework for doing that. 

Gaussian approximation. In most of the paper 'P(R) 
is chosen Gaussian of zero mean. This is an important 
restriction because the density of eigenvalues is not uni- 
versal, i.e., it depends on the specific form of V(R). Even 
worse, for complex financial systems this JPDF is typi- 
cally far from Gaussian (cf. App. for an introduction 
to non-Gaussian effects on financial markets). 

The Gaussian approximation is dictated by analytical 
difficulties — it leads to complicated solutions yet much 
more tractable than in a non-Gaussian case — hence, it 
is a natural place to begin the analysis. Furthermore, it 
does provide a basis for non-Gaussian extensions, e.g., 
to the Student t-distribution through random volatil- 
ity models, in which one considers the variances to be 
random variables (cf. App. IA 21 for an introduction and 
Sec. lIIIB~7l for an application), or to the free Levy distri- 
bution through free probability calculus (cf. Sees. Ill C 31 
IIII B 81) . 

Covariance functions. A real multivariate Gaussian 
distribution of zero mean is fully described by a two-point 
covariance function, 
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where the brackets stand for the averaging over "P(R). 
The attention in the literature has been mainly on the 
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equal-time (b = a) covariances, described by an N x N 
real symmetric matrix Cjj = (Ri a Rj a ), where stationar- 
ity (independence of a) has been assumed. This pa- 
per revolves around a more involved object, the time- 
lagged (t = b — a 0) covariances, described by an 
N x N real time-dependent matrix Cij(t) = (Ri a Rj. a +t) , 
obeying C(i) = C(— i) T , where translational symmetry 
in time (over the considered period T) has again been 
assumed. 



2. Measurement noise 

Unbiased estimators. A basic way to experimentally 
measure the above covariance matrices for a given com- 
plex system is to take a historical time series of some 
length T of the values i?, a , compute the realized covari- 
ances Ri a Rj a or Ri a Rj^ a+ t, and average them over the 
past time moments a — this results in unbiased estima- 
tors which in the matrix notation read 



1 rp 

c = — RR , 

T 



(2a) 



c(t) = — RD(«)R T , D ab (t) = 5 a+t , b . (2b) 

Weighted estimators and the EWMA. Other choices 
are also in use, e.g., "weighted estimators" in which the 
past realized returns Ri a arc multiplied by some real and 
positive weights w a , (i) the smaller the older the data 
(i.e., this sequence increases with a) which reflects the 
fact that older measurements are less relevant for today 
estimation of the covariance, and (ii) obeying the sum 
rule ^ 2o =1 w\ = l to make this case consistent with 
[([2afl. (|2bl) ]. Denoting W a b = w a S a b, this prescription is 
equivalent to R — >• RW, i.e., 
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In yet other words, this is equivalent to the standard 
estimators [ (|2a[) . (|2bj) ] with a modified true covariance, 

(RiaRjb) — > WaWb(RiaRjb) ■ 

A typical example is the "exponentially- weighted mov- 
ing average" (EWMA), 



T- 
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(4) 



with just one parameter, k £ [0,1], defining the expo- 
nential suppression of the older data with the character- 
istic time tewm a = — 1/ log ft. Estimator [ (|3a|l . Q] (with 
k = 0.94) is widespread in financial industry through the 
RiskMetrics 1994 methodology of risk assessment 

Measurement noise and thermodynamic limit. A fun- 
damental obstacle for practical usage of the estimators 
is that a description by a finite time scries will neces- 
sarily incur inaccuracies — since N time series of length 



T (i.e., NT quantities) are used to estimate the ~ N 2 
independent elements of the "true" covariance matrix, 
the error committed (the "noise-to-signal ratio" ) will be 
governed by the "rectangularity ratio" of R, r = N/T. 
In other words, c and c(t) arc consistent estimators of 
C and C(t), respectively, i.e., they tend to the true val- 
ues as r — > (i.e., for T large compared to A, which 
is a standard limit in statistics). On the other hand, a 
more practically relevant regime is when both N and T 
are large and of comparable magnitude — such as several 
hundred financial assets sampled daily over a few years — 
which is modeled by the "thermodynamic limit," being 
a standard limit in RMT, 



N,T 



r = — = finite. 
T 



(5) 



In this case, the estimators reproduce the true values 
of the covariances but with additional marring with the 
"measurement noise," the more pronounced the larger r. 

One might immediately propose to further and further 
increase the observation frequency (i.e., decrease <5i, i.e., 
increase T, i.e., decrease r) in order to suppress the noise, 
but unfortunately this procedure is limited by the Epps 
effect (cf. the end of App. IA2b|) — the very object one 
is measuring, the covariance matrix, changes with the 
observation frequency above its certain level. 

Remark that the thermodynamic limit should be com- 
bined with the following limits: 

(i) The RMT methods used in this paper are valid only 
provided that 



t«r. 



(6) 



Otherwise, finite-size effects would become relevant, as 
illustrated in Fig. [7] (a). 

(ii) For the EWMA estimators [ ([3afl . (|3b|) ]. one should 
additionally take 



1", # = T(l - k) = finite, 



(7) 



since then the characteristic time t'ewma ~ 
stretches over the scale of the entire time series. Un- 
fortunately, the limits ([6]) and ([7| arc incompatible with 
each other due to t <§; tewma, thus the theoretical re- 
sults will reproduce Monte Carlo simulations with some 
discrepancy, albeit not very big, especially well inside the 
bulk [cf. Fig. [3] (a)]. 

Noise cleaning of the equal-time covariance estima- 
tor. A basic way to observe the measurement noise is 
by comparing the eigenvalues of the equal-time covari- 
ance matrix C and its estimator c. They are impor- 
tant because they represent uncorrelated causes ("prin- 
cipal components," "explicative factors") for the corre- 
lations of the objects Ri a at the same moment of time 
(the temporal index a will thus be skipped below); in- 
deed, diagonalizing the (symmetric and positive-definite) 
matrix C, Cv^ = XjVj, yields a linear decomposition, 
Ri = Ylj=i [vj ] i e j , into uncorrelated factors of variances 
given by the eigenvalues, (ejej.) = XjSjk. This procedure 
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is known as the "principal component analysis" (PCA). A 
financial interpretation is that the eigenvectors describe 
uncorrected portfolios (investments of [vj]j part of one's 
wealth into asset i) while the corresponding eigenvalues 
quantify their risks (this is important for the Markowitz 
risk assessment method [H). 

The simplest situation is C = <j 2 \n, in which case the 
mean spectrum of the estimator c (which is known as the 
"Wishart random matrix" @) is given by the famous 
"Marcenko-Pastur (MP) distribution" Q (cf. App. ICll 
for derivation). In other words, the single eigenvalue 
A = a 2 is "smeared" by the measurement noise over the 
interval [x_,x+], where x± = ct 2 (1 ± \/r) 2 , with a den- 
sity ([C3|) . 

In the seminal papers [1, (the results below come 
from Q; cf. Fig. 1 there), the empirical spectral den- 
sity of the estimator c derived from N = 406 financial 
returns of the S&P500 index, sampled over T = 1309 
days (i.e., r « 0.31), has been fitted by the MP distri- 
bution. A fit with a 2 ss 0.74 happened to well describe 
most of the spectrum, except two regions, comprised of 
c.a. 6% of the eigenvalues: (i) One small peak around 
a very big eigenvalue Ai (c.a. 25 times larger than x+). 
(ii) Several well-separated peaks just above x+. More- 
over, the eigenvector [vi]j « 1/y/N, which represents an 
approximately equal investment in all the assets — a very 
risky portfolio, depending on the behavior of the mar- 
ket as a whole, and thus called the "market mode." In 
the first approximation, Ri sa [vij^ei, i.e., just one fac- 
tor, the market, governs the behavior of all the assets. 
Furthermore, the eigenvectors corresponding to the in- 
termediate peaks describe portfolios of strongly corre- 
lated assets from industrial sectors. To summarize, the 
bulk of the spectrum may be thought to represent a pure 
noise, whereas the spikes which leak out of it carry rele- 
vant information about the true structure of equal-time 
correlations. 

This use of the MP distribution for discerning the true 
information from noise has led to a number of "clears 
ing schemes" devised to construct a covariance matrix C 
which better reflects experimental data. For^ instance, 
in the "eigenvalue clipping" method [l(J, C consists 
of the empirical eigenvalues lying above x+ and their 
eigenvectors, while the eigenvalues below the upper MP 
edge, supposedly originating purely from the measure- 
ment noise, are all replaced by a common value A no iso 
such that the trace of the covariance matrix is unaltered. 
This scheme has been argued [ll| to slightly outperform 
classical cleaning methods such as the "shrinkage algo- 
rithm" [13 . 

These discoveries commenced a series of applications 
of RMT to econophysics. In particular, an important 
research program is to construct more realistic "priors" 
C which — unlike the trivial case above — already incor- 
porate some features of the market, then perform the 
PCA of the estimator c, and compare it with empirical 
data. More generally, one should assume a certain JPDF 
"P(R; {parameters}), dependent on parameters which de- 



scribe some aspects of the market dynamics, and compare 
the spectra of both estimators c and c(t) with empirical 
data, thereby assessing the parameters of the selected 
JPDF. In particular, quite naturally, it will be demon- 
strated that true correlations between different assets but 
at the same moment of time are more visible in the spec- 
trum of c, while c(t) becomes very useful for investigat- 
ing true time-lagged correlations. In the next Sec. II Bl a 
few such forms of the JPDF used in this paper will be 
introduced. 



B. Models of spatial and temporal correlations 

As explained above, a basic approach to distinguish- 
ing relevant information about the true correlations in 
a given complex system from the measurement noise in- 
evitably present in any covariance estimator is by com- 
paring the empirical spectrum of an estimator with its 
theoretical spectrum obtained from the assumption of 
no correlations; any difference between the two results 
then from the true correlations. This is Toy Model 1 
(Sec. IIIIBp below. A more advanced version of this pro- 
gram is to assume a certain parametric form of the true 
correlations, compute the theoretical spectrum of an es- 
timator and compare it with an experiment — this has a 
potential to more accurately determine the parameters 
of the system. This Section introduces two important 
classes of models of the true correlations which will be 
used for this purpose in the paper. 



1. Industrial sectors 

Factor models. In order to construct a prior C which 
would reflect the spike structure visible in the spectrum 
of c, a simple way is through the "factor models" ["factor 
component analysis" (FCA)]. 

For instance, a "one- factor model" ["market model," 
"capital asset pricing model" (CAPM)] (ljjj . sufficient 
to describe the market mode, assumes that each asset 
Ri (the temporal index a is omitted here) follows a cer- 
tain factor 4>q (the "market" ; a random variable with 
volatility S) with strength (3i (the "market beta"), i.e., 
Ri = f3i(f>o + Ci, where the ej (the "idiosyncratic noises," 
with volatilities o~i of comparable magnitude) are uncor- 
related with </>o and each other. The cross-covariance 
matrix reads, CV, = E 2 /3,/3j +a 2 8ij, and indeed has one 
large (~ N) eigenvalue plus a sea of (N — 1) small ones. 

More generally, a "if-factor model" [Tj-QJl further 
supposes that the idiosyncratic noises depend on K 
hidden factors, corresponding to the industrial sectors, 
&i = ^2 a =i Pia4>a + €i, where the industrial factors and 
the new idiosyncratic noises have volatilities S Q and 
Cj, respectively, and are all uncorrelated. The spec- 
trum of the cross-covariance matrix, = T, 2 /3i/3j + 
E a =i ^oft"ft«+ a i^iji ma y be shown to contain not 
only the market mode but also the K sectorial modes. 
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Other constructions along these lines exist, e.g., the 
"hierarchically- nested factor models" [13 ■ Regardless 
of the model, if the assets are Gaussian and there are 
no temporal correlations, the matrix C can simply be 
considered diagonal, with an appropriate structure of its 
eigenvalues. 

Power-law- distributed C. Another model of C has 
been suggested [llj to even better reflect the coexistence 
of small and large industrial sectors [H,[lj| — a hierarchi- 
cal, power- law distribution of the sector sizes [20| , 



Pc (A) 



fi(fi - iy (i - \ min y 



(A — 1 + /x (1 



Amin))^ 



T 0(A-A min ) (8) 



[this is the density of eigenvalues of C defined in (|B1(I ; it 
is normalized and yields -^TrC = 1], which is described 
by two parameters, a slope fi > 1 and the smallest sector 
size A m i n £ (0, 1 — l//i) [O(-) is the Heaviside theta (step) 
function] . 

In pj (cf. Fig. 2 there), N = 500 U.S. stocks has 
been observed over T = 1000 days (i.e., r = 0.5), and 
shown that the empirical spectrum of c, after removing 
the market mode, is very well fitted with the Wishart 
distribution with ([5]), where /i = 2 and A m i n = 0.35. 

This article investigates the above situation in two ver- 
sions: (i) C has K distinct eigenvalues (Toy Model 2a, 
Sec. IIII C~Tj) : (ii) C has powcr-law-distributed eigenvalues 
© with [i=2 (Toy Model 2b, Sec. IIII C 21) . 



2. Vector autoregression models 

VAR. A simple yet profound model for the evolution 
of statistical variables e.g. in macro-economy (cf. the Eu- 
ropean Central Bank's "Smets-Wouters model" (211) o r 
in finances (cf. the "linear causal influence models" J22|) 
is "vector autoregression" (VAR) in which the value of 
any variable depends linearly on the values of some other 
variables in some previous times, plus an idiosyncratic 
noise (the "error term" ) ei a which, in anticipation of the 
applications below, is chosen complex Gaussian of zero 
mean and the nonzero covariance (ei a ~e~jb) = C-f'(6 — a). 



/v 



Ri, 



(9) 



j=l b<a 



[The summation over time b is assumed to stretch from 
— oo, i.e., this is far from the beginning of the dynamics. 
The "influence kernel" Kij(c > 0) = 0.] This recurrence 
equation can be readily solved, 



N 

E E E 

n>0 di,...,a n <0 j — 1 



(10) 



[K(oi)K(a 2 )...K(a n )] i ,e 



j,a+ai + ...+a„ i 



where certain constraints must be imposed on the ker- 
nel to make this sum convergent (to be specified below 



on simplified examples). Consequently, the variables are 
complex Gaussian of zero mean and the covariance func- 
tion translationally symmetric in time (i.e., dependent on 
c = b — a) which reads 



C(c) = 



-EE- 

7l,ra>0 ai,...,a n <0 
6i,...,6 m <0 

K( ai )K(a 2 )...K(a„)- 

(m n 
b »>' - E a "' 4 
m' — l n' — l 

K(6 m )tK(6 m _ 1 )t...K(6 1 ) 



(11) 



SVAR(l). In this paper only quite drastically simpli- 
fied versions of the VAR model will be considered: (i) 
The temporal dependence in ([9]) is only on the previ- 
ous moment, K(c) = B<5 c ._i, for some N x N matrix B 
[this is abbreviated VAR(l)]. (ii) The covariance func- 
tion of the error terms is trivial in the temporal sec- 
tor, C ld (6 — a) = C ld -<5 Q f,; moreover, the spatial covari- 
ance matrix will be assumed diagonal, C- H ' = (a\ d ) 2 Sij 
[this is called "structural VAR" (SVAR)]. Under these 
circumstances, the covariance function decays exponen- 
tially with the time lag, 



C(c) 



FB^ C , 
Bl' ; lF, 



for c > 0, 
for c < 0, 



(12) 



where for short, F = J2 n > B n C id -Bt>". 

B diagonal. The matrix B will first be chosen diago- 
nal, Bij = fySij, which means that there is no correlation 
between different variables. Then Fij = fiSij, where for 
short, fi = (crj d ') 2 /(l - |A| 2 ), and there must be \/3 t \ < 1 
for all i. Hence, the covariance function becomes 



^ J (c) = / 4 |/3 l | |c| e- iar ^ c <5, 



(13) 



i.e., it is diagonal with each term consisting of a distinct 
(i) proportionality factor ff, (ii) exponential decay with 
the time lag governed by (iii) sinusoidal modula- 
tion with the time lag governed by arg(/3j) (however, the 
Pi will henceforth be always assumed real, so no oscilla- 
tions). 

Three versions of this model will be investigated in this 
article: (i) all the fii equal to each other and a\ d ' equal 
to each other (Toy Model 3, Sec. HTTP]) : (ii) all the ft 
equal to each other but arbitrary erj d - (Toy Model 4a, 
Sec IIII Elj) : (iii) arbitrary ft and a\ d - (Toy Model 4b, 
Sec. IIII E2|) . 

B with a market mode. In order to incorporate cor- 
relations between different variables, a more involved 
non-diagonal structure of B is required. For instance, 
to account f or t he market mode (cf. Sec. II B 1[) . one 
may assume (llOj that there is one asset (the "market" ) 
which depends only on itself, Ri a = ei a + ai?i ia _i, while 
any other asset depends on the market and on itself, 
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Ria = Gi a + j3R\ M ~i + jRi,a-i, f° r * > 1, where the beta 
coefficients a G (0, 1), /3 G R, 7 G (0, 1) are identical for 
all the assets. Consequently, the nonzero matrix elements 
of B read, Bu = a, Bn = /3, Bu = 7, for i > 1 (i.e., the 
diagonal and first column). The covariance function thus 
becomes 



Cii(c) = 
C l2 (c) 



/i 2 



(1 — a 2 )(Q! — 7)(1 — cry) 



/3 2 7 



1 — 7 2 \ (a — 7) (1 — 07) 

Cii(c), for c > 
Cii(c), for c < 



/3 



1 



a — 7 y 1 — a 2 
Cii(c), for c > 0, 



1 



1 — cry 



(14a) 



(14b) 



(14c) 



C u (c), forc<0 / (l-a 2 )(l-a7) 

/3 2 



C «W-(a- 7 )(l-a7) 



1-a 2 



I-7 2 



a 1+|c| , (14d) 



(14e) 



for i,j > 1 and i =/= j (i.e., five distinct matrix entries; 
the first row and column interchange upon c — >■ — c). It 
has the following (time-dependent) eigenvalues: (i) an 
(N — 2)-fold degenerate eigenvalue, 



Ai(c) 



1 - 7 2 



(15) 



describing the self-interaction of the assets; (ii) two eigen- 
values solving the quadratic equation 



II. MASTER EQUATIONS 
A. Outline 

1. Random matrix models 

The goal of this article is a basic study of the covari- 
ance estimators [([2a]). (|2b)) ] and [([5a]). ([3b)l ]. Actually, a 
few more assumptions about them will be made (with- 
out changing the notation) which do not influence the 
end results at the leading order in the thermodynamic 
limit [(0, ©], but which greatly simplify the calcula- 
tions: (i) The returns are complex instead of real, (ii) 
The indices in the Kronecker delta in the delay matrix 
(|2b[) are understood modulo T, i.e., it is unitary, (hi) The 
normalization 1 / (T—t) in (|2b[) is replaced by 1/T. Points 
(ii) and (iii) are easily justified by t <^ T; an explanation 
concerning the crucial point (i) is offered in Sec. Ill B 31 
To summarize, the following two N x N random matrix 
models are investigated: 

(i) The (Hcrmitian) "equal-time covariance estimator" 
(ETCE), 



c = -RR f . 

T 



(17) 



(ii) The (non-Hermitian) "time-lagged covariance esti- 
mator" (TLCE), 



c(t) = -RD(t)Rt, D ab (t) 



T 



>(a+t) mod T,b- 



(18) 



[An even more general model b f|29[) . which encom- 
passes the weighted TLCE (|3bj) . is also touched upon 
in Sees. lnB~2l and HTC21 ] 

In these definitions, the N xT return matrix R consists 
of complex Gaussian random numbers of zero mean and 
the covariance function of the form, 

(R ia Rjb) = Cij^b, (19a) 
(R ia R jb ) = (RiaR jb ) = 0. (19b) 



(a - 7 )((1 - a 7 ) 2 + (N — l)^ 2 ) (a 7 ) |c| - 
- ((1-7 2 ) ((a- 7 )(l-a 7 ) + (iV-l)/3 2 a)al c l + 

+ (1 - a 2 ) ((a - 7 )(1 - 07) - (N - l)f3 2 ^^ ■ 
■ (1 - «7)A 2 , 3 (c) + 

+ (1 - a 2 ) (1 - 7 2 ) (a - 7 )(1 - a 7 ) 2 A 2 ,3(c) 2 = 

(16) 

[one of them is large as 0(7V)], describing the interactions 
of the market with itself and individual assets with the 
market. This is the content of Toy Model 4c, initially 
analyzed in Sec. IIII E 31 



2. Mean spectral density 

Basic quantities. One quantity only will be calculated 
for the ETCE and TLCE — their "mean spectral den- 
sity" (MSD) [flBTJ, (|Bl8|) ]. I t is in troduced in App. [ETa] 
(Hcrmitian case) and App. IB 2 al (non-Hermitian case), 
where it is also encoded in the form of more conve- 
nient objects — the Green function [(|B3|). (|B20|) ] or the 
M-transform [([B4]). ([B23])]. 

Master equations for various true covariances. The 
rest of this Section is devoted to deriving the "master 
equations" obeyed by these objects for several cases of 
the true covariance function (|19a|) : (i) Case 1: Arbitrary 
Cij,ab (Sec. Ill B|) . (ii) Case 2: Factorized Ci^ a b = CijA^a 
(Sec. Ill Cj) . (iii) Case 3: Translationally symmetric 
in time = Cy(b - a) (Sec. Ill Dp . (iv) Case 2 + 



6 





Case 1 

[CUD, 

(fl9b|)] 


Case 2 gl]) 


Case 3 ([731 


Case 2 + 3 m 


C, A 


C = 1jv 
A = diag 


A = 1 T 


C, A 


C = liv 


Gaussian 


free Levy 


free Levy, 
/3 = 


ETCE 


[(125a). 

(l25bl 

(l27aT). 

(l27b)). (El, 

(22)] 


(1481) (recalculation) 


[([51a)l. 
JBTbl] 




[(fffa)l- 
(l78bl) (12151) 

121] 






TLCE 


[(|33a]- 
(|33hl, 
(|35al, 

(|35b I, 121, 
122|) 


'(1535}- 

(55bl). 

(35af. 

(35b)). 121, 

(221] 


((71)). bor- 
derline: 

[dZSiJ, 
d220] 


165)), bor- 
derline: 
[(pa). 

pEf] 


[(l69al). 

pb)] 


dzi 


[(f79a]-(f82l. 
121, 121] 


((95a))- 
(HSc)), 
(903), 
(90b)). 
ISB), (Pa)). 
(94b)). (E2J, 
(SI] 


[<|103a|)- 
(fTM))] 



TABLE I: Summary of the general master equations obtained in Sec. [TT] The only recalculated result is for the ETCE in Case 
2; all the rest is new. 



3: Factorized and translationally symmetric in time 
Cij, ab = CijA{b - a) (Sec. |HE|. 

All these results are summarized in Tab. |U They 
arc new discoveries except for the ETCE in Case 2 
(Sec. Ill C~T)) . which is presented for (i) completeness; (ii) 
comparison of the derivation process with the TLCE in 
the same case (Sec. Ill C~2)) ; (iii) an extension to the free 
Levy distribution of the returns (in Sec. Ill C II for the 
ETCE, then in Sec. ECSfor the TLCE). 

Application to toy models. These general equations are 
subsequently (Scc. lIII)) solved for four specific Toy Models 
(1, 2a, 2b, 3) of the true covariance function, and com- 
mented on for three other more involved Toy Models (4a, 
4b, 4c); they have already been mentioned in Sec. II Bl 
Also, it is explained how they can be a basis for exten- 
sions to non-Gaussian probability distributions of the re- 
turns: Student (Sec. IIII B 7) and free Levy (S ec, fill B 8[) . 
as well as the EWMA estimator (Sec. IIII B 6")) . 

Techniques of derivation. The method used to obtain 
the master equations is the planar diagrammatic expan- 
sion and Dyson-Schwingcr (DS) equations. All its nec- 
essary concepts are introduced in a self-contained way 
in App. IB 1 bl (Hcrmitian case) and App. IB 2 bl (non- 
Hermitian case). Occasionally (Sec. Ill C 3p . the U N- 
transform conjecture" — a mathematical hypothesis re- 
lating the eigenvalues and singular values of a non- 
Hermitian random matrix whose mean spectrum pos- 
sesses rotational symmetry around zero — is alternatively 
used which greatly simplifies calculations for the relevant 
models. 

Universality issues. Remark also that the MSD is not 
universal — it depends on the particular probability dis- 
tribution of the returns, and not only on the symmetries 
of the problem — but nonetheless it is useful in applica- 
tions e.g. to finances (cf. Sec. II A2[) . However, for the 



models with a rotationally symmetric mean spectrum, a 
universal form-factor is proposed and tested (Secs. llHB 41 
and IIII C~l")) which describes the steep decline of the MSD 
close to the borderline of its domain. 



B. Case 1: True covariances arbitrary 

In this Section, the general structure of the true covari- 
ance function [ (|19a|l . ()19bj) ] is assumed, and the technique 
of Sees. [BTb] and [Bib] applied to the ETCE and TLCE 
[actually, its generalization ((29)) ] . respectively 

1. Equal-time covariance estimator 

Step 0: Prior to that, remark that c is a product of two 
Gaussian random matrices, an d T^f (^ ne factor 

1/T could be split between the two terms in a different 
way, but with such a choice, they are of the same large- 
T order), while the method is used most easily just to 
Gaussian matrices. To remove this obstacle, a lineariza- 
tion procedure has been proposed [HI, which introduces 
the larger matrix, 

c lin ' - ( °" t ^ V (20) 
Its basic property is that upon squaring, 

(^■f=(*T A)' (21) 

it yields a block-diagonal matrix with the blocks being 
the two possible products of the two terms in question, 
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which means that c lin - has the same non-zero eigenvalues 
as c, counted twice (i.e., 2 min(iV, T)), plus \N — T\ zero 
modes. In other words, the spectrum of c lln ' consists of 
min(iV, T) 2-valucd square roots of the eigenvalues of c, 
plus \N — T\ zero modes. This translates into a close 
relationship between the holomorphic M-transforms of 
the two random matrices, 



M c n a .(z) 



2r 



1 



M c (z 2 ). 



(22) 



The cost of this linearization is the increase in the dimen- 
sion of the matrix from N x N to (N + T) x (N + T). 
The four blocks of c (and of analogous matrices hence- 
forth) will be denoted by NN, NT, TN, TT (from left to 
right and top to bottom), according to their dimensions. 
For instance, for the Green function matrix (IB2I). 



G cli „.(z) = 



qNN qnt 
QTN qTT 



(23) 



and similarly for the self-energy matrix (cf. the discussion 
of the 1LI diagrams in App. IB lb} . 

Step 1: The nonzero matrix elements of c lm ' are all 
Gaussian random numbers of zero mean, hence the full 
information about them is encoded in the propagators, 
and the nonzero ones are directly inferred from (|19ap to 
be 



1 



(24) 



Step 2: Knowing the propagators, one may explicitly 
write the second DS equation (|B13p . 



1 T 

T NN _ ±_ n y-iTT 

a, 6=1 

TT _ 1 ^ 
n ba — If 2-^1 L! i. oi,0r . 



NN 
ji ) 



(25a) 
(25b) 



and S" T = S TJV = 0. 



Step 3: Since the nonzero blocks of the self-energy 
matrix lie on its diagonal [ (|25a|) . (|25bj) ]. 



E 



TT l ' 



(26) 



the first DS equation (|B12|) implies that the Green func- 
tion matrix is also block-diagonal, with G NT = G TN = 
and 



G 



NN 



1 



qTT _ 



zl N - S 
1 



NN ■ 



zl T - S 



TT ' 



(27a) 
(27b) 



One may say that the first set of DS equations decouples 
into a spatial and temporal sector; of course, the two 



sectors are coupled by the second set of DS equations 

[CHI), <|25b])]. 

Step 4: Eqs. [(|25a|). (|25bj) . (f27ajl . (f27bj) ] are_the mas- 
ter equations for the four matrices, G 



NN 



G 



TT 



^NN 



iTT 



Their explicit solution may be attempted if the 
true covariancc function Cij^ab has been specified — six ex- 
amples will be analyzed in Sec. lIIIl Once the two blocks of 
the Green function matrix are found, adding their traces 
yields the holomorphic Green function (|B3p of c lln -, 



G c n„.(z) = 



1 



(TrG 



NN 



TrG 



TT\ 



(28) 



N + T 

which is then easily related to the holomorphic Green 
function of the estimator c [ ([22]) . (|B4|)]. and thereby to 
its MSD (lB8ll. 



2. Generalized time-lagged covariance estimator 



Consider the following generalization of the TLCE, 

(29) 



b = -RER^F, 

T 



where E and F are arbitrary constant Hcrmitian or non- 
Hcrmitian matrices (of dimensions T x T and N x N, 
respectively). [Setting E = D(i) and F = In reduces b 
to c(t).] 

Step 0: Before that, the matrix product in (|2"9")) should 
be linearized in order to obtain a Gaussian random ma- 
trix, 



^lin. 







iV 



-LR F 

VT 



Or 



(30) 



Analogously to the Hcrmitian case, the nonholomorphic 
M-transforms of the original and linearized matrices re- 
main in the simple relationship (|22p . so it is enough to 
consider only the latter. 

Step 1: The non-Hcrmitian procedure differs from its 
Hcrmitian counterpart essentially by one step, which 
is "duplication" [ (|B24[) - (|B26[) ] of the random matrix 
in question. Combining the notation from Step of 
Sec. Ill B H and (|B26j) . the Green function matrix [ (|B24[) . 
(|B25|) ] will have the block structure 



G, 



..dupiAZ.Z) = 



/qNN 


QNT 


qNN 


G NT\ 


G TN 


G TT 


G TN 


QTT 


qNN 


QNT 


qNN 


G NT 




QTT 


qTN 


QTT J 



(31) 



and similarly for the self-energy matrix. 

Step 2: The propagators (covariances) of the returns 
(|19ap allow to easily deduce the propagators of the du- 
plicated version (|B25|) of b lln - , 



^ b lin. i dupl. ] JVT [b Un.,dupl. ] TiV^ 



N T 



/v=l r=l 



(32a) 
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T T 



c=l d=l 

(32b) 



rvlm.,dupl.iWTr K lm.,dupl.iTAr\ _ 1 VVp* P Z 1 

fe=l !=1 

(32c) 



iV iV 



AT T 



[bli n. ) d U p 1 . ]r[b U,,dup 1 . ] | J V^ = l££i44 C ^c- 

(32d) 



fe=l c=l 



As explained, since the random matrix in question is 
Gaussian of zero mean, these propagators carry the entire 
information about it. 

Step 3: The second set of DS equations (|B13[) acquires 
thus the form, 

= m E EcaF k jC ik ,cbGab , (33a) 



a.b.c.k 



^ba — j, E E caFkjCik,cbG 



NN 
ji j 



(33b) 



e n_n = ]_ y E El c G tt (33c) 



a.b.c.d 



E ^4*%^™ ( 33d ) 

i,j,c,d 



S r = £ E ^Ch,»G^, (33e) 



- y E F J k F lj C kl.abG 



NN 
ji ' 



(33f) 



',b,k,i 



^ = \ E 44%^, (33h) 



i,j,k,c 



with the remaining blocks of the self-energy matrix zero. 

S^ep ^: The block structure of the self-energy matrix 
which emerges from [ (|33a|) - (|33h[) ] is such that each of 
its main four blocks (holomorphic-holomorphic, etc.) is 
block-diagonal, 



/S 7 



JJ-jlin. ,dupl 



(z,z) 



NN 





^NN 





\ 










S TT 




NN 






















/ 



(34) 



and owing to it, the first set of DS equations (|B12j) may 
be rewritten as two matrix equations (of dimensions 2Nx 



2N and 2T x 2T, respectively), 



qNN 


qNN \ 




( zl N - X NN 


_-£>NN 


qNN 


G NN J 


= 


y _^NN 


zl N - X NN 





G TT 




G TT 



zIt - S TT 









(35a) 

-l 



(35b) 



Analogously to the Hcrmitian case [([27a| ) -(|27b ]) ]. the first 
set of DS equations decouples into a spatial and temporal 
sector. 

Step 5: As in the Hermitian case, the nonholomorphic 
Green function (|B20|) of b lin - is given by (|28|l . which can 
a priori be calculated once the master equations [ (|33a|) - 
(|33TT)) ([3"5a)l (|35b|) ] have been solved. Using ||22J) and 
[ (|B23[1 . (|B21|) ] leads then to the nonholomorphic Green 
function and MSD of b. 

As explained in App. IB 2 a| another interesting part of 
the duplicated Green function matrix [ (|B24[) - (|B26[) ] is the 
order parameter — the negated product of the normalized 
traces of its off-diagonal blocks (|B3ip . 



B b un.(z,z) = - 



1 



N + T 
1 



N + T V 



TrG NN + TrG TT 
/ TrG lV W + TrG TT 



(36) 



Recall that it necessarily is a real and nonnegative num- 
ber, which zeroes only on the borderline of the mean 
spectral domain, dT> b nn. (and <92?b) (jB32[) . thus provid- 
ing a way for its analytical determination. 



3. Complex instead of real assets do not influence the 
leading order 

An argument will now be given why complex instead 
of real assets lead to much simpler models and yet do 
not change the end results at the leading order in the 
thermodynamic limit (J5j) , as anticipated in Sec. Ill ATI Of 
course, the microscopic properties for these two classes 
are altogether different, however, they are not the subject 
of this article; for an advanced calculation of some of 
such properties for a Wishart matrix with real entries, 
cf. e.g. 0, [2f|. The proof will be given only for b, since 
for the ETCE it has already been presented in [2(| , albeit 
in a less general setting of the true covariance function 
factorized (Case 2 below) — it has been demonstrated that 
all the moments of the ETCE are identical at the leading 
order for complex or real assets. 

If the random variables R ia were real, with the most 
general structure of the propagator ((T|), then in addition 
to [(pa |) -(|32d |) ]. there would be six nontrivial propaga- 
tors of b > p , namely between its blocks: NT and 
NT, TN and 77V, as well as the self-propagators of AT, 
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TN, NT, TN, e.g., 



T T 



([b li n .,dupl. ] iVT [bl in. )d up 1 . ] ^T ) = ±_Y i Y, E ~ E *£ 



ij,cd-> 



(37) 



etc. This would lead to eight additional nontrivial blocks 
in the self-energy matrix, those which are zero in Q34[) . 
e.g., 



(38) 



b,c,d,j 



etc. In other words, the solutions of all our models would 
become significantly harder. 

However, these new blocks (|38|) are one order of T 
smaller than [ (|33a| - (f33h|) ]. i.e., may be neglected in the 
thermodynamic limit. Indeed, if one follows the sum- 
mation indices in (|38p . j — > b — > c — > d, and treats it 
like matrix multiplication, one finds that is like one 
matrix entry of the resulting product, i.e., of order 1/T. 
On the other hand, in e.g. Q33a[) one has one matrix 
entry through the summation over k, plus a trace from 
a — > b — > c; this trace contributes an order of T, hence 



Y:fj N is of order 1. 



N 



For instance, in the simplest case of E = D(i), F = 1 
(the TLCE) and Ci^ a b = o~ 2 SijS a i, (no correlations; this is 
Toy Model 1 below, Sec. IIII Bp and for real assets, the 
eight entries [ (|33a[) - (|33h| ] of the self-energy matrix are 
identical as in the complex case and read 



X NT =a 2 ±(r> T G TN " T 
T \ 



T 



. 2 I(g^d) t , 



J}NT = a 2± Q TN 



T 

S =a - [G ) , 



(40d) 
(40e) 
(40f) 
(40g) 
(40h) 



but it is clear that they are one order of T smaller than 
the other entries, and thus irrelevant in the thermody- 
namic limit. 



C. Case 2: True covariances factorized 

In this Section, a special form of the true covariance 
function (jl9a| is investigated — factorized into a spatial 
and temporal part, 



C 



ij S*-ba ) 



(41) 



where C and A are some constant Hermitian matrices 
of dimensions N x N and T x T, respectively, describing 
the decoupled spatial and temporal covariances. (Note 
the convention about the indices of A.) 



X NN = [ -Tr< 



z, =a* I -Tr (D« 
£ TT = a 2 (^TrG w 



- cr | — TrG ) ljv, 



( -TrG iViv ) 1 T , 



I — Tr ^D T G TT ) ) 1 N , (39g) 



(39a) 
(39b) 
(39c) 
(39d) 
(39e) 
(39f) 



(39h) 



however, the eight other entries (|38|) are now nonzero, 

1 / ^ r TAT\ T 



1 



s™=^(G^r, 

T 



(40a) 
(40b) 
(40c) 



_Z. Equal-time covariance estimator 

Planar diagrammatics derivation. For the ETCE, the 
second set of DS equations [ (|25a|) , (|25b() ] takes on the 
form 



S TT = Ag, 



where for short, 



7^Tr(A 
.9^Tr(C 



TT\ 



NN\ 



(42a) 
(42b) 



(43a) 
(43b) 



In other words, the matrix structure of the self-energy 
blocks is discovered — they are proportional to C or A, 
respectively, with unknown proportionality constants 7, 
g, which depend on the blocks of the Green function ma- 
trix, and this by multiplying these blocks by A or C, 
respectively, and taking the traces normalized by 1/T. 

Plugging [ (l42a|) , (|42b|) ] into the first set of DS equa- 
tions [ (|27ajl . (|27bj) ]. and using the definition of the holo- 
morphic Green function matrix (|B2[) . one obtains 



qNN _ 



1 = -G C (- 

zl N - G7 7 V7 



(44a) 



10 



iTT 



1 



zIt — 



G, 



(44b) 



Now, one should combine Eqs. [ (|44a|) . (|44bp ] with the 
definitions of 7, g [ (|43a|l . (|43b[) ] in a twofold way: 

Firstly, insert the former into the latter, then use 
the definitions of the holomorphic Green function 
(|B3[) and M-transform (|B4p , which allow to rewrite 
^Tr(AG TT ) = ^Tr(A(zl T - Aff)" 1 ) = M A (z/g)/g, 
and analogously in the spatial sector. This yields two 
scalar equations for two scalar unknowns 7, g, 



[ - \ = rM c ( ^ 



(45) 



Secondly, take the traces of [ (|44ap . (|44bl) ] and substi- 
tute into [(]2"8"]). (J22)], using also This implies that 
the holomorphic M-transform of the ETCE is expressed 
through 7, g as 



1 



M c (z 2 ) = - 73 



(46) 



The master equations [ (j45j) . (|46|) ] may be rewritten 
in an even more concise way by exploiting the notion of 
the "holomorphic TV-transform" which for any Hermitian 
random matrix H is the functional inverse of the holo- 
morphic Af-transform, 

M H (JVh(*0) = A H (M h (z)) = z. (47) 

This quickly leads to one simple equation for M = M c {z) 1 

rMN A (rM)N c (M) = z. (48) 

It is the only master equation presented in this paper 
which has already been known (cf. e.g. [2H 127|). 

Free probability derivation. Equation (148p can also be 
derived in another even simpler way [281 ] . based on free 
probability theory of Voiculescu and Speicher [29], l3Cj . 
Even though this derivation is already known, it will 
now be sketched because of three reasons: (i) It is based 
on an application of the multiplication law (|49jl to a 
product of two Hermitian matrices, a procedure which 
may a priori not be valid. For the ETCE, it does work. 
However, an analogous reasoning may be applied to the 
TLCE as well and it is quite surprising to discover that it 
yields a wrong result, albeit only slightly (cf. Sec. Ill C"3j) . 
(ii) It can easily be adjusted to the situation of the re- 
turns distributed according to the free Levy law, which 
is a model of a non-Gaussian behavior of financial assets 
(cf. App. IA 1 dp . However, this extension has been ac- 
complished [31[ only for C = ljv, A = It and the free 
Levy distribution with zero skewness. Therefore, it will 
in what follows be worked out for arbitrary C, A and 
the Levy parameters a, /3, 7 (except the case a = 1 with 
/3 7^ 0). (iii) This result will in turn constitute a basis for 
a derivation of the master equation for the TLCE and 
free Levy assets (cf. Sec. IIIC3]) . 



In order to rederive (]4"8"]) . recall the "multiplication 
law" of free probability. If Ui and U2 are two free uni- 
tary random matrices ("freeness" is a generalization of 
the notion of statistical independence to noncommuting 
objects; here it may be thought of simply as indepen- 
dence), then the MSD of their product is obtained by 
multiplying their A-transforms, 



N- 



UiU 2 



CO 



z + 1 



N Vl (z)N V2 (z). 



(49) 



This formula may be generalized to matrices other than 
unitary, although caution must then be exercised. For 
instance, it may be applied to Hermitian matrices, even 
though their product is gencrically non-Hermitian. It 
does works for certain models such as in [H, HH or the 
below derivation, however, Sec. Ill C 31 contains a coun- 
terexample. (Cf. the end of App. I A fdl for the free prob- 
ability addition law in the Hermitian world, and the end 
of App. IB 2 al in the non-Hermitian world.) 

Now if the true covariance function of the returns is 
factorized (j4"Tj) . they can be recast as Re C 1 ^ 2 RiA 1 / 2 , 
where the new returns arc IID of variance 1. Combining 
cyclic shifts of terms with the multiplication law ([49]) 
allows then to relate the ETCE c for the original returns 
R to the ETCE C! for the IID returns Ri, 



N c (z) = 



A iT3 , 4R t c (z) 



cycl. 

4- 

""iRiARjC^ 
{49) 

cycl. 

I rz* 



(l + z)(l + rz) 



-N 



^rJri {rz)NA(rz)Nc(z) 



(50a) 



cycl. o 

i_ rz z 

(l + z)(l + rz) 



N Cl (z)N A (rz)N c (z). (50b) 



Finally, inserting N Cl (z) = (1 + z)(l + rz)/z, which 
readily follows from the MP holomorphic Green function 
(IC2p . reproduces the master equation (|48| . 

Generalization to free Levy returns. In the free Levy 
instead of Gaussian world, formula (|50bp remains intact 
[actually, (|50a|) will be used; also, the normalization is 
rather T~ 2 / a , cf. (|A5|) ] . only a new expression for N Cl (z) 
needs to be worked out. 

It can be accomplished by the "projector trick" [28] ], 
which is to define a rectangular (N x T; assume here 
N < T) free Levy random matrix as the projec- 
tion with P = diag(ljv, Ot-n) of a square (T x 
T) free Levy Hermitian matrix L. This implies 
T _2//q R{R 1 = LPL, in which changing cyclicly the or- 
der of terms and using the multiplication law (|49p (again 
for a pair of Hermitian matrices; it works this time as 
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well) yields N T _ Va - R i R (z) = (z/{z+l))N P (z)N Ij 2(z) 
= ((z + r)/(z + l))Nv{z). 

Secondly, simple algebra applied to the definition 
(1B231) yields M L2 (z 2 ) = \{M^(z) + M L (-^)). Now, for 
z in the upper half-plane, Ml(z) follows from (|A8|) . while 
Ml(— 2) from the same equation only with j3 — > —(3. 

Combining these results leads to the following set 
of two master equations for two complex unknowns, 
M = M c (z) and auxiliary 5, 



(rM + 1 + S) a (rM + 1 - 5)°> 



LP ' rM + S 
r 2 M 2 ^1 - 
1 

' r 2 M 2 _ 5 2 C 

where for short, 

,4 



(rM + l) 2 

L7ra _ 

2 ' 



rM -5 

N A (rM)N c (M) 



7 



e™ ap , for a G (0,1), 

e i7r(a-2)^ j for a g ^ lj2 ). 



(51a) 



(51b) 



(52) 



In particular, for a symmetric distribution ((3 = 0, i.e., 
ip = 1), Eq. (|51ap is trivially satisfied with 6 = 0, upon 
which Eq. (|51b() turns into 



j™(rM) 2 ( a -V N A (rM) a N c (M) a 



r 



(53) 



For a = 2 and 7 = 1, it further reduces to the Gaussian 
result 



2. Generalized time-lagged covariance estimator 



For b (|29[) . one is able to follow the steps in the above 
planar diagrammatics derivation to a certain extent, rec- 
ognizing similar structures, however, a reduction to a sin- 
gle equation like (|48|) is yet to be accomplished. 

Really, the only manageable simplification so far is to 
write the second set of DS equations [ (|33a|) - (|33h[) ] in a 
form analogous to [ (|42a) . (|42bjl . (|43a| . (|43bj) ]. 



(54a) 
(54b) 



where 



7i 






^2 


72 y 




5i 






hi 


92 J 





bTr 



( CF 7l 


Cm \ 


\ FtCF77 2 


Ftc 72 J ' 


( AE S i 


A/11 \ 


y E^AE/lg 


EtA ff2 J ' 


AEG TT 


EtAEG TT N 


AG TT 


EtAG TT , 


CF QiViV 


FtCFG^ 17 


CG NN 


FtCG™ 



(55a) 



(55b) 



These along with the first set of DS equations [ (135a[) . 
(|35b|) ] form a set of master equations in the considered 
case. 



In order to compare them with their ETCE counter- 
parts, define 



E up 




E 



down 



It 








Et 



(56) 



and recall A duph = diag(A, A) (|B25j), p lus analogously 
in the spatial sector. Then, [ (|35a[) . ()35b[) ] combined with 
[ (|Ma|) . (|54b)) ] read 



1JV _ 



G T = 



qNN 


qNN 


qNN 


qNN 


G TT 


g tt\ 


G TT 


qTT J 



1 



jndown(^dupl.pjpup ' 

(57a) 

1 



Z — E down A du P L CJE u P : 



while the definitions [ (|55a|) . (|55b[) ] become 

; bTr (A du P L E up G T E down ) 



G = 



7i 




h 


V2 


72 J 




9i 






h 2 


92 j 





-idupl.-pup^iV-pdownA 



(57b) 



(58a) 



(58b) 



One recognizes now similarities between [ (|54a|) . (|54b[) ] 
and [ (|i2a|) . (l42bl) ]. [ (IBTal) . (l57bl) ] and [ (pa) , (jiib]) ]. 
[ (I58al) . Ipbll ] and [(pajl. (|4lh|)]. 

A way to reduce these master equations to a simple 
single equation like (|4"5| for the ETCE has unfortunately 
not been found. The reason may essentially be traced 
back to the fact that the block trace (|B27j> — unlike the 
usual trace — does not exhibit the cyclic property; this 
obscures a way to write an analog of Eq. (|46|) . which re- 
lates the spectral information about the estimator to the 
constant matrices T, Q. In other words, just as (|48| uti- 
lizes the notion of the holomorphic iV-transform (|4"7| (of 
the Hcrmitian matrices C and A), one would expect the 
desired simplification to be possible only once a proper 
definition of the iV-transform for non-Hcrmitian matri- 
ces has been found, which is not yet the case (cf. [34[ for 
recent progress in this direction). 



3. Rotational symmetry 

Rotational symmetry of the MSD. There is an im- 
portant property, often encountered throughout this pa- 
per, which the MSD of a non-Hcrmitian random matrix 
X (e.g., the TLCE) may possess — rotational symmetry 
around zero, i.e., dependence only on 

R = \z\, (59) 

in which case it is sufficient to consider just its radial 
part, 

p% d -(R) = 2wR 0x(*,Z)| w=a = 

_d_ (60) 

dR 



— M x (z,z)\ M=R , 
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where the lower line originates from [ (|B23[) . (|B21[) ] . 

N -transform conjecture. In such a situation, a very 
helpful tool is the following hypothesis which has been 
put forth and extensively tested in [HI, [H, Hf| HH : 

(i) For any non-Hermitian random matrix model X 
whose nonholomorphic M-transform (and therefore, the 
MSD) is rotationaUy symmetric around zero, 

M x (2,z)=OTx(i? 2 ), (61) 
one may invert (|6ip functionally, 

9Kx (ttxC*)) = z, %c(Wt x (i? 2 )) =R 2 , (62) 

which is named the "rotationaUy symmetric non- 
holomorphic A-transform." 

(ii) For the Hermitian random matrix X^X, the holo- 
morphic A-transform is the functional inverse of the 
holomorphic M-transform (|47p. 

(iii) The TV-transform conjecture claims that these two 
quantities arc related through 

iVxtxW = ^x(4 (63) 

It is useful because often one of the models, X or X^X, 
is more tractable than the other (concerning computation 
of the MSD), and therefore Eq. ([63|) gains access to that 
more complicated model. 

Relationship between the TLCE and ETCE assuming 
the rotational symmetry. A general procedure to discern 
whether the MSD of the TLCE for a given true covari- 
ance function has the rotational symmetry will not be 
attempted in this work. Potential implications of this 
property will nonetheless be analyzed and two specific 
examples outlined. 

Application of the above hypothesis to the TLCE re- 
quires caution concerning the order of terms: 

Firstly, introduce the random matrix 
c(t) = iR*RD(t), which differs from c(t) only 
by the order of terms and dimension; conse- 
quently, their nonholomorphic M-transforms relate 
as M e(t) (z, z) = rM c(t) (z, z) = rM. 

Secondly, assume the rotational symmetry around 
zero (to be verified a posteriori), and write the rota- 
tionaUy symmetric non- holomorphic iV-transform (|62[) . 
m m {M m (z,z)) = R 2 . 

Thirdly, consider the Hermitian random ma- 
trix c(i)1c(i) = D(t)tc 2 D(t), where c = ^RtR. 
The A-transform conjecture (|63|) reads 

Ar c(t)tc(t)( 2; ) = (( Z + 1 )/ Z )^c(t)( Z )' 0n the 0ther 

hand, changing the order of the constituents and using 
the unitarity of D(i), one finds M 6 ( t )tc(t) ( z ) = Mc 2 i z )- 

Fourthly, recall M-^z 2 ) = ±(M 6 (z) + M- c {-z)). 

Fifthly, changing again the order of terms, 
M e (z) = rM c (z), where c is the ETCE (fTT|) . 

Collecting all the above formulae leads to an equation 
which relates the M-transforms of the TLCE and ETCE, 



M = ReM c ±iRj- 1 + -T7 • (64) 



[The expression under the square root is always nonneg- 
ative, cf. the discussion below. The sign, reflecting the 
two possible square roots, is irrelevant.] 

Single ring conjecture. Another hypothesis for non- 
Hermitian random matrices with the rotational symme- 
try (|6ip has been proposed [35l |36| — which may be re- 
garded as an extension of the "Feinberg-Zee single ring 
theorem" [3Tl - f40| — that their mean spectral domain is 
either a ring or an annulus; the radii of the enclosing 
circles are denoted by i? e xt. and i?; n t.- Assuming this is 
true, for the TLCE, it is a general caveat that M = 
on the external circle and M = — 1/r on the internal 
one, where the latter exists only for r > 1. Indeed, 
on the borderline, the nonholomorphic and holomorphic 
M-transforms coincide. But the only holomorphic func- 
tion compatible with the rotational symmetry (i.e., de- 
pending on R) is a constant, M c m (z) = Mq. More- 
over, it is known from RMT that in the external out- 
side (i.e., containing z = oo), the holomorphic Green 
function must behave as G c u)(z) = (Mo + l)/z ~ 1/z, 
for z — > oo — which implies Mq = 0. On the other hand, 
in the internal outside (i.e., containing z = 0), the value 
of M is related to the zero modes: The MSD (|B2ip 
reads, ^c\G c ( t ) (z) = (M + l)S (2 \z,z), where the rep- 
resentation of the complex Dirac delta (|B19p has been 
used. Now, zero modes appear in the spectrum of the 
TLCE only if N > T, and their number is (A - T), 
i.e., their density, (1 — l/r)S^ 2 ' (z,z); indeed, the ma- 
trices RD(i)R^ and R^RD(t) have the same nonzero 
eigenvalues, while the larger one has additionally | — T| 
zero modes. This implies Mq = —1/r in the internal out- 
side, which exists only for r > 1. To summarize, M is 
a smooth function of R which monotonically grows from 
max(— 1/r, —1) at R = R m t. to at R = R cx t.- 

Example 1: Trivial temporal covariances. It will be 
proven (cf. Sec. [HE3| that the MSD of the TLCE cer- 
tainly exhibits the rotational symmetry in the case of 
the true covariance function factorized (|4ip with arbi- 
trary C and A = It- The master equation for the 
holomorphic M-transform of the ETCE (|48|) becomes 
then (rM c (z) + l)A c (M c (z)) = z, therefore Eq. $64$ 
turns into a complex equation for two real unknowns, 
M = M c ( t }(z,z) and auxiliary m, 

M + im = M c [ ±R\ - { 1 + — ] ^ r I • 

y y V rM J rm - i(l + rM) I 

(65) 

Remark that if C has a finite number of distinct eigen- 
values, this can be recast as a set of polynomial equa- 
tions (cf. Toy Models 1 and 2a below). For C with an 
infinite number of distinct eigenvalues, this is typically 
more complicated (cf. Toy Model 2b). 

The radii of the circles enclosing the mean spectral 
domain are obtained from (|65p by expanding it around, 
respectively, M = or if r > 1 also M = —1/r, which 
both are consistent with assuming an expansion around 
m = 0. The external one is given by the first and sec- 
ond moments of C (|B5b|l . while the internal one by its 
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holomorphic M- and TV-transforms, 



i c 1 + rm c ,2, 



^nt. = rNci-l/rfM'c (N c (-l/r)) Q(r 



(66a) 
(66b) 



[Remark that (|66bp may be recast in a shorter but per- 
haps less useful way, R? nt = —2r/d z Nc(z)~ 2 \ z =-\/ r .] 

Example 1 — a surprisingly wrong derivation using the 
multiplication law. Before proceeding to a second ex- 
ample of a rotationally-symmetric TLCE, recall that the 
master equation (|48[) has been rederived (cf. the mid- 
dle part of Sec. Ill C lj) using the free probability multi- 
plication law @5J) by relating (|50b|> the ETCE for true 
covariances factorized into C, A (|4"Tj) to the ETCE for 
IID returns. One may therefore ask whether the master 
equation (|65|) could be obtained in a similar fashion by 
relating the TLCE with the spatial covariance matrix C 
to the TLCE c x (t) for IID returns (which is Toy Model 
1 below). The derivation should start from R = C 1//2 Ri 
and proceed as 



cycl. 

i 



f63l 



z + 

cycl. 



J N C Cl {t)f Cl (t)c( z ) 



z + 



——N C 2(z)N Cl(t)tci(t) (z) 



|63j 

4- z 



jN C 2(z)m cl{t) (z) 

(generally incorrect!). 



(67a) 



Now the rotationally-symmetric nonholomorphic N- 
transform of Ci (t) is given by functionally inverting (|120|) . 
which leads to the following equation 



M = M C 2 R 



1 



rM 



rM(l + r + 2rM) 2 
(generally incorrect!). 



(68) 



It is surprising to discover that it is different from 
Eq. ([6"5)) . It seems that the only place the above deriva- 
tion may have failed is the application of the multipli- 
cation law to a product of two Hermitian matrices, as 
explained in Sec. Ill CH It is even more surprising to find 
out that the right and wrong master equations produce 
solutions with only a very slight difference in the MSD; 
Figure [7] (b) confirms the agreement of Monte Carlo sim- 
ulations with Eq. (I65[) and their disagreement, albeit very 
small, with Eq. (1551) . 

Example 1 — generalization to free Levy returns. An- 
other point is to recall that the master equation (|48|) for 
the ETCE for factorized true covariances (|41[) has been 



generalized (cf. the last part of Sec. Ill C 1[) to the situa- 
tion of free Levy returns [ (|51a|) . (|51b|) ]. These formulae 
[with A = It, i.e., Ma(z) = l/(z — 1)] can be combined 
with ([64]) to yield a set of two master equations for the 
nonholomorphic M-transform of the TLCE for arbitrary 
C in the free Levy case, 



<P- 



(rM + 1 + 5 + irm) a (rM + 1-6 + irm) c 



rM + 5 + irra rM — 5 + irm 

((rM + 1 + irm) 2 - 5 2 ) (M + im) 



(69a) 



(rM + 1 + irm) 
1 

r 2 (M + im) 2 -5 2 ' 



1 



R a 

7 



N C (M + im) 



(69b) 



where real m and complex S are auxiliary unknowns. 

In particular, for zero skewness, Eq. (|69a|) is trivially 
satisfied with S — 0, while Eq. (|69b[) becomes 



(rM + 1 + irm) (M + im) 



(M + im) 



. e i7r f r a ~ 2 



N C (M + im) 



(70) 



R" 

T 2 



Further setting a = 2 and 7 = 1 brings it back to the 
Gaussian result (f6"5j) . 

A derivation of the values of the external and internal 
radii of the borderline of the mean spectral domain will 
not be presented in full generality, just commented on for 
special values of the parameters in Sec. MI B 81 

Example 2: Trivial spatial covariances and diagonal 
temporal covariances. Another situation when the rota- 
tional symmetry is apparent from an inspection of Monte 
Carlo data is the factorized covariance function (I4T1) with 



C = ljv and arbitrary diagonal A = diag(crf, . . . , ct t ). 
The master equations for the TLCE [<[5ia |l -((5BE| ) . (|55ajl . 
(|35b[) ] happen to be very cumbersome [especially the ma- 
trix inversion in (|35bj) ]. therefore a proof will not be 
offered, rather the assumption will be made and nu- 
merically verified a posteriori. Since Eq. (|48|) becomes 
r(M c (z) + l)N A (rM c (z)) = z, Eq. (J64J) turns again into 
a complex equation for two real unknowns, M and aux- 
iliary m, 



M+im = -M A l±R 




rM J r(m-i(l + M)) 
(71) 

Recall that this equation holds conjecturally for A diag- 
onal, not arbitrary. 

The radii of the enclosing circles are found in an anal- 
ogous way as above to be 



R,, 



rm A,i + r2 ™A,2, 



f? 2 — 



(r-lf 



+ ( r - l)WA,-2 



(72a) 

6(r-l), (72b) 
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where the negative moments are defined analogously to 
()B5b[) but with negative powers, and remark a different 
placement of the r factors in (|72aj) as compared to ()66al) . 

One could also directly generalize Eq. (|7Tj) to the free 
Levy case; the result is not needed in this article and will 
not be presented. 



Remark finally that in a situation when the entries 
B a b themselves depend on T, the above procedure may 
not work. This is the case e.g. for the EWMA estima- 
tors [([3a)l. (|3b]). gl)] in the limit ©, to which a different 
method needs to be applied, namely Eq. (JTTJ) based on 
the rotational symmetry of the MSD (cf. Sec. IIII B~6l) . 



D. Case 3: True covariances translationally 
symmetric in time 

Another interesting special form of the true covariance 
function (|19a|) is when its dependence on time is transla- 
tionally symmetric, 



2. Equal-time covariance estimator 



Cij : ab = Cij(b — a). 



(73) 



Note that the right-hand side may be treated as a time- 
dependent N x N matrix, C(c); it obeys C(c) = C(— c)L 



1. Fourier transformation 

Since the thermodynamic limit ([5|) is assumed, it 
proves very convenient — especially with the symmetry 
([73")) present in the system — to encode all the T xT ma- 
trices (T — > oo; assume that the temporal matrix in- 
dices range from — oo to +oo) by their Fourier transforms. 
Namely, instead of any such matrix -B a b, use 



-t-oo 



B(p,q)= J2 ^ i(ap ' bq) Bab, 



(74) 



— oo 

where p, q S [— 1/2, 1/2], which conversely reads, 

-1/2 r-1/2 



B„ 



dpdqc 



-2iri(ap—bq) 



1/2 J -1/2 



B(p,q). (75) 



Remark that in particular, the Kronecker delta S a b is 
mapped to the Dirac delta S(p — q), while the matrix 
multiplication B = B1B2 translates into the integration 
B(p,?) = X! / 1 2 /2 dr5i(p,r) J B 2 (r,g). 

If a matrix satisfies B^ = B(a — b), where B is some 
function of time, then its Fourier transform is propor- 
tional to the Dirac delta, B(p 7 q) = 8{p — q)B(p), with 



+00 



B(p) = J2 ^ 1CP B(c), 

C— — OO 

,1/2 

B{c) = / dpe- 2nicp B{p). 
J-1/2 



(76a) 
(76b) 



In particular, matrix multiplication is simply the 
multiplication of these reduced Fourier transforms, 
B(p) = Bi(p)B2(p). This language of Fourier transforms 
for T x T matrices will henceforth always be used pro- 
vided that the temporal translational symmetry (|73j) 
holds. 



The general master equations for the ETCE [ (I25al) . 
()25bj) . (|27ap . (|27bjl ] will now be written in the situation 
([73f in the Fourier language. 

The second set of DS equations acquires the form, 



1/2 



= / dpC(p)G TT (p), (77a) 
J-1/2 

£ TT ( P ) = ±Ti-(c(p)G NN ), (77b) 

while in the first set of DS equations, the spatial sec- 
tor is unchanged and the temporal sector is Fourier- 
transformed, becoming a scalar equation, 



qNN _ 



1 



G TT (P) 



zl N - £ 
1 



NN ' 



Z 



E TT (p) 



(78a) 
(78b) 



Finding a solution to [ (|77aj) . (|77b) . (|75a|) . (|78b[) ] may be 
attempted once C(c) [and thus C(p) (|76a|) : it is Hcrmi- 
tian for any p] is known. 



3. Time-lagged covariance estimator 

Consider now the TLCE [not the more general ran- 
dom matrix b (|2"5|) . since arbitrary E, F may break the 
translational symmetry in time]. 

The second set of DS equations [ (|33a|) - (|33h[) ] can thus 
be rewritten as follows: In the spatial sector, 



,NN 



■\NN 



f 1/2 dpe- 2 ^C(p)G TT ( P ), (79a) 
J-1/2 



/ 
1/2 

-1/2 
1/2 

-1/2 
1/2 



d P C(p)G J 1 (p), 



d P C(p)G J 1 (p), 



dpe 2 * itp C(p)G TT (p), 



(79b) 
(79c) 
(79d) 



'-1/2 

while in the temporal sector, 



E-(p) = -Tr(C( P )G 



NN 



(80b) 
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NTT 



(p) = -Tr(C(p)G 



.VA" 



(P) 



e 2^itp Tr 



C(p)G 



WW 



(80c) 
(80d) 



In the first set of DS equations [ (|35a|) . (|35b|) ]. the spa- 
tial sector is unchanged, while the temporal sector is 
Fourier-transformed, becoming a 2 x 2 matrix equation, 
featuring explicitly doable 2x2 matrix inversion, 



qNN 


qNN 


qNN 


qNN 



( ,i s^NN 




\ 


_J}NN 


ZLjsr — 2j 


G TT (P) 


G TT (P) \ 




G TT (P) 


G TT ( P ) ) ' 





(81a) 



z-± TT ( P ) 


-t TT {p) 


-s tt (p) 


z-t TT {p) 





' z-t TT (j>) 


£ tt (p) 


W T (p) 1 




z - t TT {p) 



(81b) 



where for short, 
W T (p) = 



NTT 



(p) 



(p)) 



(82) 



>:'"( P )>:"(7;). 



E. Time-lagged covariance estimator for Case 2 + 
3: True covariances factorized and translationally 
symmetric in time 

Recall that in the case of factorized true covariances 
(|41[) . the master equations for the ETCE have been sim- 
plified to a single equation (|48p (cf. Sec. Ill CTj) , which has 
not been accomplished for the TLCE (cf. Sec. Ill C 2|1 . 
Some simplification becomes however possible if in ad- 
dition the temporal translational symmetry (|T3|) is as- 
sumed, 



Cij,ab = CijA(b - a). 
1. Arbitrary C and A 



(83) 



Firstly, the temporal sector of the second set of DS 
equations [ (|5Da|) - (|80dl l] acquires the form 



NTT 



NTT 



NTT 



(p) 



tp A(p)h, 



( P ) = A(p)t 2 , 
{ P ) = A{p)t 3 , 
(p) = e 2 " H *A(p)U, 



(84a) 
(84b) 
(84c) 
(84d) 



where the information about the model is now encoded 
in four complex parameters, 



h = ^Tr ( 
h = ^Tr ( 
h = ^Tr ( 



NN\ 



S 4 = ^Tr ( CG" J 



(85a) 
(85b) 
(85c) 
(85d) 



Notice that the definition of the Green function matrix 
[ |B24l . (|B25]) ] implies i 4 = h; it will moreover be shown 
[cf. Eq. ([TPD]) ] that 



h = -t 2 t 3 G 



h>0. 



(86) 



where h = only on the borderline of the mean spectral 
domain (|10ip . Hence, there really is one complex and 
one real unknown. 

Secondly, substitute [ ()84a|) - (|84d[) ] into the temporal 
sector of the first set of DS equations (|81b[) . 



G TT ( P ) 


G TT ( P ) 


G TT ( P ) 


G TT (P) 



W T (p) 



z- e 27Titp A(p)U 


Mp)h 


A{p)t 3 


z - e- 27ritp A(p)ti 



(87) 



where (|52"]l. 
W T (p) = 



A( P ) ( e - 27ritp ^i 



ztA 



+ A{pf (tit 4 - t 2 t 3 ) 



Thirdly, knowing (the Fourier transforms of) the tem- 
poral blocks of the Green function matrix [ ([87|) . (|88|) ] , one 
may compute the spatial blocks of the self-energy matrix 

(89a) 
(89b) 
(89c) 
(89d) 
where for short, 

f 1 ' 2 1 . 





= C (zh- 






= Ct 2 I 2 , 






= Ct 3 I 2 , 




jyNN 


= C (zh+ 





(90a) 
(90b) 



(Note, I\ + = while I 2 is real.) Remark that the 
matrix structure of the spatial blocks of the self-energy 
matrix is discovered — they are all proportional to C. 
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Fourthly, thanks to this proportionality, calculating 
the spatial blocks of the Green function matrix (|81a|) only 
requires inverting a 2 x 2 matrix, 



G 
G 
G 
G 

where for short, 



NN _ 


zl N - 


C {zh+ 


- Hh) 






W N 




NN _ 


Ct 2 I 2 








W N ' 






NN _ 


Ct 3 I 2 








W N ' 






Wn _ 


zl N - 


C (zh- 


- Uh) 



w 



(91a) 
(91b) 
(91c) 
(91d) 



W 



jv 



= \z\ 2 - 



C (I 2 (zh + zU) - z 2 I 1+ - z 2 h-) + 
C 2 f/| (hU - t 2 t 3 ) - I 2 (ztj^ + ztth-) + 



(92) 



(Here all the matrices commute, so they may be treated 
as numbers.) 

Finally, in order to obtain a set of equations for the 
basic unknowns, £1,2,3,4, the results [ (|91a|) - (|91dj) ] should 
be plugged back into the definitions [ (|85ap - (|85dj) ]. which 
yields 



h 
h 

where for short, 



zJi — (zli- 

t 2 I 2 J 2l 
t 3 I 2 J 2l 

z J i — (~zli- 



- hh) J 2 , 



- Uh) J 2 



J i = -Tr 



I, 
T 
1 



W N ' 

c 2 



Tr- 
T W N 



(93a) 
(93b) 
(93c) 
(93d) 



(94a) 
(94b) 



Note that I\±, I 2 , Ji, J 2 all depend — generically in a 
very complicated way — on the unknowns £1,2,3,4- 

Nonholomorphic solution. As a general caveat 
(cf. App. |B2|) . the master equations [ ([93a|) - ()93dp ] will 
have a nonholomorphic solution — valid inside the mean 
spectral domain, z € T> c ( t yin. , and describing the MSD 
there — and a holomorphic one — for the outside of the do- 
main, whose matching with the nonholomorphic solution 
provides the equation of the borderline of the domain. 

In the nonholomorphic sector, there must necessarily 
be t 2 ,3 ^ [cf. the discussion around Eqs. (fTUII)) . ([TUT]) ]. 
Simplifying [ (|93ap - (|93d|) ] accordingly yields 



J 1 = 



Zlo 



■h = — , 



(95a) 
(95b) 



z 2 h+ = z 2 h. 



(95c) 



Once £1,2,3,4 are found from [ (|95ap - (|95cp ] — which, as 
mentioned, seems to be a formidable task — the non- 
holomorphic Green function of the linearized estimator 
c(t) hn - fl3n|) is given by an analog of (|2"5|) . with the re- 
spective traces calculated from (|91a|) and (|87p . 



(zh+ - txl 2 ) Ji, 



-TrG TT =zI -Uh + , 



where for short, 



T - *T 1 



1/2 x 
-1/2 



(96a) 
(96b) 

(97a) 
(97b) 



These quantities can be expressed through J±, J 2 
and Ii±, I 2 , respectively, by using yTrl^v = r and 
f-l/ 2 dp = 1, 



2 Jo 



= r- {h {zh +zt 4 ) - z 2 I 



1+ 



Ji- 



\z\ 2 I = l + zhh- + Ztih+ - (hU ~ ht 3 )I 2 . 



IaihU-tata) 
+ \z\ 2 I 1+ I 1 _)j 2 , 



(98a) 
(98b) 



Inserting [ (|98a|) . (|98bp ] into [ (|96a|), (|96b|) ] . and then into 
the analog of (|28j) . using also (|B23p . and finally into the 
analog of (|22p — we express the desired nonholomorphic 
M-transform of the TLCE (in the argument z 2 ) through 



^1,2,3,4, 

M c(t) (z 2 ,z 2 ) = -(zhh- - (ht 4 - t 2 h)h)- (99) 

Eqs. [H25U-G33, pa) . (l90bj) . ([88]), pijl . (IMbl) . (pjl : 
([9^)l ] thus conclude the derivation of the master equations 
under the considered circumstances (|83|) . 

Holomorphic solution. One may also evaluate the or- 
der parameter [ (|B3ip . (|36|) ] inside the mean spectral do- 
main by using the identities (|87p (its off-diagonal equa- 
tions) and [ (|9lb)) . (|9Tc)) ]. 



B c{t) nn.(z,z) = h ^-^ 



I2J1 + I1 



where for short, I\ = J_^, 2 dpA(p)/W T (p) 



(100) 

_ x /2 ^xjj^yjj,! r, Since the 

left-hand side is real and nonnegative, and so is the 
bracket on the right-hand side — so must h be (|86p . 
Therefore, the equation of the borderline f|B32|) reads 



0. 



(101) 



Special cases. In the remaining part of this Section, 
two subcases of ([83| will be investigated: cither C = 1 jv 
(Sec. lnE~2ll or A(b - a) = 5 ab (Sec. IITE31) . 
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2. Trivial C and arbitrary A 



As a first subcase, suppose 



L/y. 



A(b — a) = arbitrary. 



(102) 



Then, the spatial traces [ (|94a|) . (|94b[) ] greatly sim- 
plify, and so do the master equations [ (|95a|) - (|95c|) ]. be- 
coming Ii + = ~z/z, Ii_ = z/~z and I<i = r / (tit± — £2^3)- 
Moreover, (|99p turns into (here it proves simpler to use 
the nonholomorphic Green function instead of the M- 
transform) G c ( t )(z 2 , z 2 ) = t\/(rz). Collecting all these 
results (and changing the argument z 2 to z) yields the 
final set of master equations, 



V2 d — -L_ 

-1/2 P Hp) ~ \G\ 2 + h' 



1 



1/2 1 

dp — — 

1/2 w(p) r 2 A(p)e 2mt P 



where for short, 



w(p) 



G - 



.,2711/ /I 



rA{p) 



+ h, 



(103a) 
(103b) 

(104) 



where the unknowns are: complex G = G c i t ){z, ~z) 
and real and nonncgativc h = h c (t)(z, z), with 
h c (t){z 2 , z 2 ) = h/(r 2 \z\ 2 ). Notice that there generi- 
cally is no rotational symmetry here. 



3. Arbitrary C and trivial A 

As a second subcase, suppose on the contrary, 

C = arbitrary, A(b — a) = 6 a b. (105) 

Recall that the master equation in this situation has al- 
ready been derived (|65|) using the A^-transform conjec- 
ture; now it will be rigorously proven. 

Nonholomorphic master equations. Immediately, 
A(p) = 1, and the integrals [ (|90a[) . (|90b[) ] can be per- 
formed explicitly: Introducing a new variable 



c- 27dtp , i.e., 



1/2 1 

M---) = — 

-1/2 



1 



27ri -/C(0,1) 



d«-(...) 



(106) 

[notice that the dependence on t is lost — this subcase is 
independent of the value of the time lag, of course as 
long as t <C T, cf. Fig. [7] (a)], and applying the method 
of residues — the integration requires finding roots of a 
quadratic polynomial in v, 



V(v) = ~zt lV 2 + (\z\ 2 + t x t A + h) 
which leads to 

II .1 



zU, 



(107) 



h+ = 



2ni 



C(0,1) 



vV(v) 



1 flz^ + hU + h 



h- 



2zt± 
1 

27Ti 



dv 



C(0,1) 



V(v) 



1 flz^ + hU + h 



2zh 
1 
2?ri 



-1 , 



dv 



C(0,1) 



V(v) 



(108a) 

(108b) 
(108c) 



where for short, 



s = \f (\z\ 2 + + hf - 4\z\ 2 t!t4. (109) 



[The two roots of V(v) (|107p obey |«i||«2| = 1, hence, one 
of them lies inside the centered unit circle C(0,1), and 
one outside it; only the interior one contributes to the 
integrals.] 

Consequently, Eq. (|95c|) implies ti = f/z, t4 = f/~z, 
with / e M. The other two master equations [ (|95a|) . 
()95b[) ]. turn finally into two real equations for two real 
unknowns f,h, 



1 



-Tr- 
T W 



s, 



(110a) 
(110b) 



where substituting the results [ (fT08a|) - (fT08c)l ] into ([92]) 
yields 

W N = B.+C(C+2f)^ (R+ f2 ~ R *~ hR ) - (HI) 
and recall that in this notation Eq. (|109|) reads 



R+ J —+h) -Af 2 



(112) 



Once these equations [ (|110a|l - (|112|) ] are solved (after 
specifying the true spatial covariancc matrix C), then the 
nonholomorphic M-transform (|99|) is expressed through 
/, h as 



M = M c(t) {z,z) 



1 R 



2r 



f 2 



- 1 



(113) 



In [ (|110a|) - (lll3[) ] the argument has already been 
changed from z 2 to z. Moreover, all these identities de- 
pend only on the radius R = \z\, so the rotational sym- 
metry has been demonstrated. Conversely, one may say 
that deviations from rotational symmetry in the MSD of 
the TLCE point toward nontrivial true temporal covari- 
ances. 

Nonholomorphic master equations — a simpler form. 
The above master equations can be considerably sim- 
plified. One should start from rewriting J\ and J2 [i.e., 
the left-hand sides of pTUajl . (|110b[) with (flTT]) ] in the 
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language of the holomorphic M-transform of C, which is 
done by factorizing W N = w(C - Q(C - (), 



Ji 



J 2 = 



Mc(0-M c (C) 

^(c-c) ' 

^M c (() -(M C (Q 



(114a) 
(114b) 



This leads to 



M + im = M, 



c 



./ 2 - (h + s)R-R' 



2VkR 
/ 2 + (fe - a)fl + i? 2 
2ri?/s ' 



, (115a) 
(115b) 



which is supplemented by [ (|112|) . f|l 13[) ] . Using these last 
three equations, one expresses /, h. s through M, to, 
and substitutes them into (|115ap , which finally yields a 
concise complex equation (|65|) for two real unknowns, M 
and auxiliary to. 

Borderline of the mean spectral domain. It will now be 
proven that I? is a ring for r < 1 and an annulus for r > 1 , 
with the radii [ (|66a|l . (|66b[) ]. To find the borderline, the 
original form of the master equations is more convenient, 
in which one should set h = (|10ip . To begin with, 
Eq. (|112j) becomes s= \R— f 2 /R\, which inserted into 
Eq. (|113p confirms the two possible values that the non- 
holomorphic M-transform of the TLCE may acquire on 
the borderline, M = or M = — 1/r. Furthermore, ma- 
nipulating Eqs. [ (fTTTJaj) . (|110b|) . ([TTT]) ] with h = leads 
to two solutions: 

(i) External circle: If |/| < R (which corresponds to 
M = 0), the master equations become 



/ = rmc.i, 
R 2 -f 2 = rm c , 2 , 



(116a) 
(116b) 



which reproduces (|66a|) . 

(ii) Internal circle: If |/| > R (which corresponds to 
M = —1/r; there must be r > 1), the master equations 



become 



R 2 i T 
T = r iv Tr - 



/ 



i 



f-R 2 = ri-Tr- 
N 



which, using the definitions [ (|B3 
(l66bl) [with / = -JV c (-l/r)]. 



(117a) 
(117b) 

B4|) . (gll)], reproduces 



III. TOY MODELS 



Outline 



Toy Models 1, 2a, 2b, 3 



In this Section, the general master equations devel- 
oped in Sec. [TT] (and summarized in Tab. |T| — which are 
a gateway to the MSD [ ijEIjl . (|Bl8|) ] of the ETCE ^ 
and TLCE (fl~8|) for arbitrary complex Gaussian assets — 
are applied to three specific models of the true covari- 
ance function (|19aP which are designed to approximate 
certain aspects of the real-world financial markets and 
other complex multivariate systems, according to the in- 
troduction in Sec. II Bt 



(i) Toy Model 1 (Sec. UlTB]) : The assets are IID. This 
is the simplest case; nevertheless, it is useful as the "null 
hypothesis" (cf. Sec. II A 2]) . i.e., any discrepancy between 
the spectra derived from experimental data and from this 
model signifies genuine correlations present in the system. 

The MSD of the TLCE has first been calculated 
in [4l], |42j, but their formula is demonstrated to be 
wrong (Sec. MI B 5p . and the right one is presented 
(Sees. |niB21IITriT3liniB4j ): it is noteworthy that the 
only analytical result to date concerning the spectrum 
of the TLCE was incorrect — it reveals how involved the 
subject is. 

This result is generalized to the EWMA estimator for 
Gaussian assets (Sec. MI B 6p . as well as the standard 
estimator for two models of non-Gaussian behavior of 
the returns: Student (Sec. MI B 71 it comes in two ver- 
sions: with a common random volatility and IID tempo- 
ral random volatilities) and free Levy (Sec. MI B 8P dis- 
tributions. 



(ii) Toy Model 2 (Sec. InTC]) : The assets have arbi- 
trary variances but no temporal correlations; they come 
in two versions (cf. Sec. II B ip : either a finite number 
of variance sectors (Toy Model 2a, Sec. MI C ip . or vari- 
ances distributed according to a power law with a lower 
cutoff (Toy Model 2b, Sec. IIIIC 2p . This and the fol- 
lowing models signify another than the null hypothesis 
approach which is to build more and more realistic toy 
models, comparing them to data, thereby verifying their 
validity and assessing their parameters. 



(iii) Toy Model 3 (Sec. IIII D[) : The assets are indepen- 
dent from each other but exhibit exponentially-decaying 
autocorrelations (identical for all the variables) , which is 
a simplest version of the SVAR(l) model (cf. Sec. II B 21) . 

It is enlightening to compare the results for the TLCE 
and ETCE, hence the MSD of the latter is rederived and 
illustrated in App.lClfor Gaussian assets and Toy Models 
1 (App.[CT|), 2a (App.EH and 3 (App.lCl. 
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TABLE II: Summary of the main results obtained in Sec. Mil and App.lClfor the ETCE and TLCE: the pertinent (i) Section, 
(ii) Figures, (iii) formula for the MSD or master equations from which it derives, (iv) equation of the borderline of the mean 
spectral domain (for the TLCE). 



2. Toy Models 4a, 4b, 4c 

Two reflections naturally arise after analyzing the 
above models — firstly, they are still quite far away from 
financial reality; secondly, the analytical derivations of 
the MSD are already increasingly complicated. Conse- 
quently, it is a valid research program to both construct 
more realistic toy models and to enhance the techniques 
presented in this paper (especially on the numerics side) 
to be capable of handling such models. Even though the 
second part is beyond the scope of this article, three more 
models are announced and initially analyzed (i.e., the 
pertinent master equations are explicitly written down 
or commented on, and the Monte Carlo simulations are 
presented; Sec. MI E)) — in contrast to the above examples, 
they are mixes of both the spatial and temporal correla- 
tions from Toy Models 2a and 3, being more advanced 
examples of the SVAR(l) class: 

(i) Toy Model 4a (Sec. IIIIE1|) : There is a number of 
variance sectors and exponentially-decaying autocorrela- 
tions identical for all the variables. 

(ii) Toy Model 4b (Sec. IIIIE2) ): There is a number of 
variance sectors and exponentially-decaying autocorrela- 
tions different for each sector. 

(iii) Toy Model 4c (Sec. IIIIE3j) : A SVAR(l) model 



with a simple nondiagonal covariance function originat- 
ing from the market mode. 

For these examples also the ETCE is analyzed because 
for models 4b and 4c the true spatial and temporal corre- 
lations are not factorized, unlike in the remaining cases, 
which causes the need to use the new master equations 
developed in Sec. Ill Bl 

The main results obtained for the above six models are 
summarized in Tab. [Til 



B. Toy Model 1 (null hypothesis): Uncorrelated 
Gaussian/Student/free Levy assets and 
standard/EWM A estimator 

1. Definition 

The simplest instance of the true covariance function 
(|19ap . serving also as the null hypothesis, is the factorized 
version (|4T]) with 

C = a 2 l N , A = 1 T , (118) 

i.e., the variables Ri a are IID Gaussian random numbers 
of zero mean and some variance a 2 . 
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FIG. 1: Monte Carlo eigenvalues versus theory for the TLCE for Toy Model 1. 



2. Borderline of the mean spectral domain 



(One easily checks that R[ n t. < Rcxt. for all r. 



For the TLCE, the considered structure of the true co- 
variance function falls under the discussion in Sees. Ill LTsI 
and lll C"3l The mean spectral domain is either a disk (for 
r < 1) or an annulus (for r > 1), with the radii of the 
enclosing circles [ (|66a|) . (|66b[) ]. 



3. Mean spectral density 



R 2 ^=a 4 r(l+r), 



lit, 



r - 1 



-6(r- 1). 



(119a) 
(119b) 



The relevant master equation is (|65p with 
Mc{z) = o~ 2 /{z — a 2 ), and elimination of m leads 
to a cubic polynomial equation for the nonholomorphic 
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FIG. 2: Demonstration that the result for the MSD for the TLCE for Toy Model 1 based on the Abel transformation is incorrect. 



M-transform, M = M c ^(z, z), 

4r 3 M 3 +4r 2 (l + 2r)M 2 + 

+ r ((l + r)(l + 5r) - R 2 a ) M+ (120) 

+ (l+r) (r(l + r) - i? 2 ) = 0, 

where for short. R a = R/a 2 , and R £ [Rmt., -Rext.]- The 
rotational symmetry around zero is evident. Once 
M is found from (Tl2"0"l) , the (rescaled) radial MSD 
p™ d - ee a 2 ff c ^j(R) = dM/dRv (JHOJ). Actually, this rela- 
tionship allows to write a depressed cubic equation di- 



rectly for the MSD, 

= a(/ CT ad -) 3 + Vr i -+c, (121a) 
a = Ar 3 ^- (ll + 14r + 2r 2 ) i? 2 - 

-(l-r 2 )(l-r) 2 ), (121b) 

b = r( - Ri +2(1 + r)(5 + r)i? 2 - 

-(1-r 2 ) 2 ), (121c) 
cee 2(1 + r) 2 R a . (121d) 
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FIG. 3: Monte Carlo eigenvalues versus theory for the EWMA-TLCE for Toy Model 1. 




FIG. 4: Monte Carlo eigenvalues versus theory for the TLCE for Toy Model 1 and two versions of the Student t-distribution 
of the assets. 



One may check that the discriminant of (|121a|) is always 
negative, i.e., it has one real root which is the desired 
MSD, 




(122) 



where for short, b = b/a and c = c/a, and two redun- 
dant complex conjugate roots. Recall also that this case 
is independent of the time lag t [cf. (j 106(1 ] . This is 
the time-lagged counterpart of the MP distribution (|C3[) 
(cf. App.IUTI). 

Figures [T] [(a), (d)] depict, respectively for r = 0.1 
or r = 10, the numerically-generated eigenvalues of the 
TLCE and the theoretical borderline(s) [ pT^ - p!9b)) ]. 
which are seen to accurately enclose the mean spectrum, 
with an exception of a small number of eigenvalues leak- 
ing out (this being a finite-size effect). In Figs.[T][(b), (e)], 
formula (|122p is shown to precisely reproduce the his- 
tograms of (the absolute values of) these numerical eigen- 
values, respectively for r = 0.1, 0.5, 0.9 or r — 10, 2, 10/9, 
again except the vicinity of the borderline(s), to which a 



finite-size reasoning needs to be applied. 



4- Universal erfc scaling 



As explained in App. IB 1 b[ the method of planar 
diagrammatic expansion employed in this paper works 
for infinite matrix dimensions (JSJ). However, for any 
non-Hcrmitian random matrix model X (of dimensions 
N x N) whose MSD exhibits rotational symmetry around 
zero (cf. Sec. Ill C 3|) there has been put forward the "erfc 



conjecture" 3j|33j,l35j,l36| which claims that the univer- 
sal finite-size behavior of the MSD close to the borderline 
of the domain (i.e., as it seems, always a ring or an annu- 
lus) is obtained simply by multiplying the density with 
the following form-factor, 



N,q b 



Rb;Sb( R ) = ^ CrfC (<?b S b (R - Rb) 



(123) 



for each circle (one or two) R = i?b which constitutes 
the borderline. Moreover, the sign Sb is +1 for the 
external borderline and —1 for the internal one; qt, is 
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FIG. 5: Monte Carlo eigenvalues versus theory for the TLCE for Toy Model 1 and the free Levy distribution of the assets. 



a parameter (one for each circle) depending on a par- 
ticular model, whose derivation will not be attempted 
(this requires genuinely finite-size techniques), but its 
value adjusted by fitting to Monte Carlo data; finally, 



crfc(x) 



^= dy cxp(— y 2 ) is the complementary er- 
ror function. 

This universal form-factor has been first calculated for 
the Ginibre unitary ensemble (|B35|) [IHIIJ] (cf. also [isf). 
then for a product of two rectangular Gaussian ran- 
dom matrices with IID entries ITuI : in both cases, the 
mean spectral domain is a disk, and the parameter q has 
been analytically determined. Inspired by this perfor- 
mance, the result has been conjectured and tested nu- 
merically for a product of an arbitrary number of rectan- 
gular Gaussian matrices [32|,[33| (one borderline), as well 
as a weighted sum of unitary random matrices from the 
Dyson circular unitary ensemble (CUE) [H| (one or two 
borderlines), or an arbitrary product of these last two 
models J3(| (one or two borderlines). 

In Fig. Q] [(c), (f)], numerical tests of this hypothesis 
are presented, with excellent agreement; in particular, 
the steep decline of the MSD close to the borderline(s) is 
reproduced. 



5. Previous solution based on the Abel transformation is 
incorrect 



The MSD of the TLCE in the situation (fTTgjl has first 
been evaluated in [U |42j . It will now be demonstrated 
that their solution is incorrect. 

Abel transformation. The main idea of fill |42| 
is conceptually similar to the A-transform conjecture 
(cf. Sec. Ill C 3[) — which recall constitutes one way to de- 
rive the above formula (|122p in which one relates the 
rotationally symmetric MSD of a non-Hcrmitian ran- 
dom matrix X to its "absolute value squared" X^X. 
On the other hand, the authors of [4l|, [42j conjecture 
that under the same assumption of the rotational sym- 
metry, one may relate the spectra of X and its "real 



part" (X + X''')/2 through the so-called "Abel transfor- 
mation," which in both ways reads 

P(x+xt)/ 2 (V2x) =- [°°dR ^ {R \ , (124a) 



Pi (R) = - 2-R / dx 



y/R 2 - X 2 

akP(X+Xt)/2 (V2x) 



Vx 2 - R 2 



(124b) 



(generally incorrect!). 



They then reinforce this hypothesis by verifying that 
[ (I124al) - (1124b[) ] work for X being the Ginibre unitary en- 
semble with some variance a (|B35|) . in which case the 
real part is the Gaussian unitary ensemble (|B9j) with 
variance c 2 /2 (this may be easily checked by computing 
the pertinent propagator). Plugging both MSD [ (|B47[) . 
(|B17|) ] into the above formulae, one proves that they are 
indeed related by the Abel transformation. 

Arguments against [ \124aty -(124b)]. This coincidence 
is however accidental. The authors of [I?} have first fal- 
sified the Abel transformation conjecture by constructing 
a counterexample, X = H1H2, where Hi^ are two inde- 
pendent GUE random matrices. 

This conjecture does not work for Toy Model 1 ei- 
ther. To show it, notice two facts: Firstly, formula (|122|) 
is certainly correct, as it has been obtained not only 
through the A-transform conjecture but also the dia- 
grammatic expansion method (which is a solid proof), 
plus it perfectly agrees in the bulk with the Monte Carlo 
data (cf. Fig. [T]). Secondly, the MSD of the real part 
of the TLCE is known (cf. e.g. [28| for derivation) — the 
holomorphic Af-transform M = -M( c (t)+c(t)t)/2( z ) obeys 
a quartic polynomial equation, 

r 2 M 4 + 2r(l + r)A/ 3 + (l + 4r + r 2 - z 2 ) M 2 + 



2 1 



— \M 



1 = 0. 



(125) 



Substituting these results into [ (|124a|l - (|124bD ] . one finds 
a discrepancy between the left- and right-hand sides of 
both equations, as demonstrated in Fig. [21 
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To summarize, it is fortunate that the MSD of this 
model is given by a solution (|122|) to a simple depressed 
cubic equation, and not by a much more complicated 
combination of the quartic equation (|125p and the inverse 
Abel transformation (|124b[) . 



6. Generalization to the EWMA estimator 

A first extension to be considered of the TLCE in the 
situation (fTTg)) (let for simplicity a = 1) is its EWMA 
version [(|3b]l. in the limit (0. As noted in Sec. II A 21 
this limit requires t ~ T, which however is incompati- 
ble with the method used — therefore, the results will not 
reproduce Monte Carlo simulations precisely. 

The master equation is the conjectured (|71| with di- 
agonal A = W, since as explained in Sec. Ill C 31 the 
exact master equations seem hard to solve. The limit 
([7| allows to pass in it from a discrete temporal descrip- 
tion (summation over a) to a continuous one (integra- 
tion over a/T G [0, 1]) since w\ — > 
which yields an explicit expression for the holomorphic 
M-transform of A, 

M A (z) = 

= I(log((l-e*)z + tf)- ( i 26) 
-log((l-e*)* + 0e*)). 

It is straightforward to numerically construct the MSD 
from this master equation. 

The borderline is a ring (r < 1) or an annulus (r > 1) 
with the radii [ ([72a)) . ([72b])]. 



7. Generalization to the Student t- distribution 
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They tend to the standard values [ f|119a|) . (|1 19b)) ] when 
0. 

Figures[3][(a), (b)] present the MSD derived from [ ([71)1 , 
(|126p ] and compare it with Monte Carlo simulations, for 
r = 0.5, three values of if = 2,4,6 and, respectively, 
t = 1 or t = 100. The concord is satisfactory for t = 100, 
since this value is large enough to be comparable with 
tewma and small enough for the finite-size effects not to 
be manifested. On the other hand, for t = 1 there are 
anticipated discrepancies, especially close to the border 
of the domain. They may be checked to decrease with 
decreasing $ (if $ is less than about 1, one may quite 
safely use small t). Moreover, one may verify that for 
Tewma T and t <C T the agreement between theory 
and numerics is excellent (not shown here). One may also 
infer from these plots that for decreasing d they approach 
the MSD of the standard TLCE, while larger $ means a 
shorter effective length of the time series, hence a more 
"smeared" density. 



As explained in App. IA 1 cl a simplest model of a non- 
Gaussian behavior of financial assets assumes that they 
are uncorrelated and distributed according to the Stu- 
dent t PDF ([A3)) . Moreover, it is known (cf. App. lA2a)) 
that a Student random variable may be thought of as 
a Gaussian random variable whose volatility <j is itself 
a random variable such that 1/cr has the gamma dis- 
tribution ()Allj) . This prescription will now be applied 
to transform the MSD of the TLCE to the Student case 
from two versions of the Gaussian model, both with inde- 
pendent assets and either (i) a common random volatility 
a for all the assets [HI, |49[ ; or (ii) IID random volatilities 
a „ at distinct time moments, but identical for all assets 
i gj. 

Version 1: Common random volatility. Assume first 
that the Ri a arc IID Gaussian with volatility er, as above 
([Tig]) Multiplying the MSD (fT2"2")) with the weight func- 
tion (|A1 1|) and integrating it over a (and changing the 
integration variable to £ = R a ) leads to the following ra- 
dial MSD of the TLCE, 



StucL.rad. / t->\ 
Pc(t) ( R ) 



0i' 



2»/ 2 r (f ) i?!+^/2 



(128) 



for R > (i.e., it is unbounded for all values of the 
parameters), where the limits of integration are [ (|119al) . 
(fll9b]t ]. 

This function is plotted and tested against Monte 
Carlo simulations in Fig. [D (a), for r = 0.5, 9 = JJi and 
three values of \x = 3, 5, 10, with excellent concord in the 
whole domain — since there are no boundaries, the erfc 
form-factor is unnecessary. Even though the spectrum 
spreads to infinity, it decays fast and is effectively more 
restrained to the region around zero than for Gaussian 
assets (formally corresponding to fx oo), the more so 
the smaller [i is. 

Version 2: IID temporal random volatilities. A more 
realistic approach would however be to choose the i?, a to 
be independent Gaussian with arbitrary volatilities o~ a , 
dependent on a but not i (jAlOj) . and eventually to assume 
these to be IID random variables distributed according 
to ()Alip . In this case, analogously as for the EWMA 
estimator (cf. Sec. IIII B 6)1 . the master equation is the 
conjectured ([7T)) with 



M A (z) 
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(129) 



and it can be straightforwardly solved numerically — with 
proper care, as (|129[) cannot be expressed by elementary 
functions— to produce the MSD of the TLCE. 

The borderline is a ring (r < 1) or an annulus (r > 1), 
extending to infinity for /i < 4, with the radii [ (I72al) . 
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They tend to the standard values [ (|119a|) . (|119bp ] in the 
Gaussian limit, 9 = yjjl, /i — > oo. 

The solutions have been found and checked with Monte 
Carlo data in Fig. 0] (b), for r = 0.5, 9 = ^[JL and three 
values of /i = 3, 5, 10, with perfect agreement, except of a 
small area close to the border (for fj, = 5, 10, i.e., when the 
MSD is bounded), to amend which the erfc form- factor 
should be used (not shown here). Note by comparing 
Figs, (a) and (b) that although the qualitative behavior 
of the MSD is analogous for Versions 1 and 2, its precise 
form is quite different. 



8. Generalization to the free Levy distribution 



(yet not big, especially for no skewness). It is also notable 
that numerical histograms are nearly identical for zero 
and nonzero skewness, while the theoretical results differ, 
especially close to zero. Moreover, one may argue (by 
expanding M and m around zero, cf. Sec. Ill C 3[) that 
the domain in these cases extends to infinity, so only a 
part of it is shown. 



C. Toy Model 2: Gaussian assets with distinct 
variances 

In this Section, effects of the presence of nontrivial 
true spatial covariances will be analyzed; since there are 
no temporal correlations here, it is enough to consider a 
diagonal prior C. 

Only the standard TLCE and Gaussian assets will be 
investigated below, albeit generalizations to the EWMA 
version of the TLCE or to Student or free Levy assets — 
as for Toy Model 1 (cf. Sees. IIIIB6( IIHB7I HUB 8ft— are 
within reach of the methods described in Sec. Ill C 31 



Another non-Gaussian model considered in this paper 
is the free Levy law (cf. App. IA 1 djl . The master equa- 
tions for the TLCE in the situation of the true covari- 
ance function factorized (|41[) into arbitrary C and trivial 
A = It has been developed in Sec. Ill C 31 to be [ (I69al) . 
()69b|) ] or (|70|) for zero skewness. Setting in them C = ljy 
(i.e., JVc(z) = 1 + l/z; let a = 1 since the variance is 
generalized by 7) simplifies them to 

(rM + 1 + S + irm) a _ (rM + 1-5 + irm) a 



rM + 5 + irm 



rM -5 



(131a) 
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with real m and complex 6 auxiliary unknowns, or in the 
case of zero skewness, 



(rM + 1 + irm)(M + 1 + im) 
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(132) 
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Figure [5] shows Monte Carlo histograms and the MSD 
following from numerical solutions to these equations for 
r = 0.5 and the Levy parameters a = 3/2 (cf. the end of 
App. IA 1 d[) . 7 = 1 and respectively (3 = (a) or f) = 0.8 
(b). The numerical eigenvalues have been obtained from 
Wigncr-Levy random matrices, as being more relevant for 
finances, hence some expected discrepancies are observed 



1. Toy Model 2a: Variance sectors 

Definition. To begin with, in the fashion of the factor 
models (cf. the first part of Sec. II B lj) . consider a finite 
number K of distinct variances, 



C = diag( a\ 



A = In 



Ni times A/2 times Nk times 

(133) 

where denote pk = Nk/N, for k = 1,2,..., K, which 
obey J2k=iPk = 1) an d which are also assumed finite 
in the thermodynamic limit. The holomorphic M- 
transform of C, which is the basic ingredient of the mas- 
ter equation, thus reads 



A" 



M c (z) = £ 



Pk 



(134) 



fc=i K 



Borderline of the mean spectral domain. As explained 
in Sees. Ill C 31 and HIE 31 the MSD is rotationally sym- 
metric around zero, and its domain is either a disk (r < 1) 
or an annulus (r > 1). The external radius (|66a[) reads 
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The internal radius (|66b[) is a function of a solution / 
of a polynomial equation of order K, 
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FIG. 6: Monte Carlo eigenvalues versus theory for the TLCE for Toy Model 2a. 



This equation has no positive solutions (/ > 0, re- 
quired for positive Rf nt ) for r < 1, and exactly one 
positive solution 2 for r > 1. Indeed, the function 
•F(f) = X)fc=i J0p — 1/r has the first derivative always 
negative, and singularities (vertical asymptotes) at the 
points — in other words, J-(f) is piecewise decreas- 
ing in the intervals (— oo, — cr^}, [—<^K'~ a K-i\i ■■•> 
[— af,+oo). A first implication is that all its zeros are 
real. Moreover, T(0) = 1 — 1/r, which is negative or zero 
for r < 1 and positive for r > 1; and J-"(+oo) = —1/r, 
which is always negative. This completes the proof. 

Mean spectral density. The master equation is (|65|) 
with (|134p . and it is straightforward to solve it numeri- 



cally for any K. 

Figurcs[n][(a), (b), (c)] show comparisons of the numer- 
ically solved master equation with Monte Carlo data, all 
perfectly supporting the theoretical results. Figures (a) 
and (b) differ by having r = 0.1 and r = 10, respectively, 
but both assume K = 2 with a\ = 1, erf = 5 and three 
values of p\ = 0.1, 0.5, 0.9. In Fig. (a), one may recognize 
two "lumps," corresponding to the two variances smeared 
with the measurement noise, of sizes proportional to the 
relative multiplicities. However, this is rather an excep- 
tion due to a small value of r and a considerable dif- 
ference between the variances; in general, unfortunately, 
any structure originating from separate variances will be 
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FIG. 7: Two tests for the Monte Carlo eigenvalues and theory for the TLCE for Toy Model 2a: the spectrum is not rotationally- 
symmetric for large t (a); the derivation of the master equation based on free probability is surprisingly incorrect (b). 



lost under the noise, cf. Fig. (c), where three values of 
r = 0.02,0.1,0.5 are considered for K = 5 with of = 1, 
of = 2, of = 4, erf = 6, erf = 8 and equal relative multi- 
plicities 0.2. In other words, these plots reveal that the 
MSD of the TLCE is far less sensitive to true spatial co- 
variances than it is the case for the ETCE (cf. App. IC 21 
and Fig. 125))— it is thus the ETCE which should be the 
main probe of true spatial covarianccs. 

Universal erfc scaling. Since the MSD is rotationally 
symmetric around zero, the erfc hypothesis (| 123[) can 
be applied, finding accurate numerical confirmation de- 
picted in Figs. [6] [(d), (e), (f)]. 

Miscellaneous. Two facts announced above can now be 
verified numerically. First, it has been said that t needs 
to be much smaller than T for the planar diagrammatic 
method used in this paper to work (cf. Sec II A 2)1 . Indeed, 
Fig- 0(a) shows the Monte Carlo eigenvalues of the TLCE 
for Toy Model 2a, with K = 5, of = 1, of = 2, erf = 4, 
a\ = 6, <t| = 8, equal relative multiplicities 0.2, as well 
as r = 0.5 or r = 2 (in the inset), and t = 250, i.e., of 
magnitude comparable with T = 1000. In both cases, 
the rotational symmetry of the MSD is lost, even though 
it has been proven for t <C T in Sec. HIE 31 Actually, 
one may infer from these (and other not shown) plots 
a conjecture that the mean spectral domain for i ~ T 
has rather a (T/t)-fold symmetry around zero, polygon- 
shaped for r < 1 and star-shaped for r > 1. 

Second, in Sec. Ill C 31 an alternative derivation of the 
master equation (|65[) has been presented — based on the 
free probability multiplication law (149[) applied to a prod- 
uct of Hermitian matrices, which a priori may not work, 
and inspired by an analogous derivation for the ETCE 
(cf. Sec. mnH)— which resulted in Eq. ([68]). Figure[7](b) 
compares, for the same values of the parameters as in 
Fig. (a) except for t = 1, a Monte Carlo histogram with 
numerical solutions of these two equations, demonstrat- 
ing that Eq. ((65)) works perfectly, while Eq. ((68)) leads 
to some discrepancy, albeit very small (it seems to be 
smaller for decreasing r). It deserves a deeper explana- 
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FIG. 8: Monte Carlo eigenvalues versus theory for the TLCE 
for Toy Model 2b. 



tion why the multiplication law produces a wrong result 
and why these two quite different equations lead to very 
close solutions. 



2. Toy Model 2b: Power-law-distributed variances 



Definition. Suppose now that the variances are dis- 
tributed according to the power law {8} with /i = 2 
(cf. the second part of Sec. II B 1[) . The holomorphic M- 
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transform of C is straightforward to calculate, 
1 



M c (z) = 



(1 - 2A min + zf 

■ ( (1 - 2A mi „ + z) (- (1 - 2A min ) 2 + z) + 
+ 2(1- A min ) 2 z- 

■ ( log (z - A min ) - log (1 - A min ) - nij^ . 

(137) 

Borderline of the mean spectral domain. As usual with 
the rotational symmetry, the domain is a disk for r < 1 
or an annulus for r > 1. The external radius (|66aP for an 
arbitrary value of /i reads 

2 J 00 > for M < 2, 

oxt ' ~ \ r (r + 1 + (1 - A min ) 2 j^j , for n > 2, 

(138) 

i.e., for /i < 2 (which encompasses the model considered 
here) the eigenvalues of the TLCE spread to complex 
infinity [cf. (pOa|) ]. 

For the internal radius (|66b[) (/.t = 2 and r > 1 as- 
sumed), the following expression is obtained, 



p2 
-"-int. 



2/ 2 (/-/+)(/-/-) 
(/ - /o) (/ + A min ) ' 



(139) 



where for short, 

V - 1± 



± V(r - 1) (r - 1 + 8A min (1 - 2A min ))) , (140a) 
/o = 1 - 2A min , (140b) 

and where real / > — A m ; n is a solution to the nonlinear 
equation 



o = Hf) 



= 2(1- A min ) log - 

/-/o 



(141) 



[1 — 2A m j n ) 



r/ 

■(-(r-l)(l 
-(2 + r-4A min )/ + / 2 ). 

Observe now the following properties of J^if)'- (i) 
It decreases inside the interval [/_,/+] and increases 
outside of it (notice that f±, /o are real and 
-A min < /_ < 0, /+ > f > 0). (ii) It has infinite 
limits, lim^_ A + J-(f) = — oo, limy_j. - J~(f) = — oo, 
lim/_ K) + J"(/) = +oo, lim/_> +00 F(f) = +co. (iii) 
■F(fo) = 0. All this implies that F(f) increases from — oo 
at / = —A m in ( a vertical asymptote) to some value (a 



maximum) at / = /_, then decreases to — oo at / = 
(a vertical asymptote); on the other side of this asymp- 
tote, it decreases from +oo to some value (a minimum) at 
/ = /+, then increases to +00 as / — > +00. Since /o > 
is a root of this function (but uninteresting as it yields 
an infinite R? nt ) , it means there is exactly one other root 
fi > 0; it leads to a positive R? Dt , and therefore is the 
right solution. Moreover, there would be two roots in the 
interval (— A mm , 0) if and only if there were ^F(f-) > 0; 
however, this may be verified to never be true. To sum- 
marize, Eq. (|141|) has always precisely one real and posi- 
tive solution fi besides /o which yields a real and positive 

Mean spectral density. The master equation is (|65p 
with (I137p . and with proper care can be solved numeri- 
cally. 

Figure [8] compares its only solution yielding a positive- 
definite and normalized MSD to Monte Carlo simulations 
for A m i n = 0.35 (cf. the end of Sec. II B lj) and three values 
of r = 0.1,0.5,0.9. It is notable that the agreement is 
perfect in the whole domain — there is no upper edge, 
thus no need to use the crfc form-factor. 



D. Toy Model 3: Gaussian assets with temporal 
covariances 

1. Definition 

In the last example, there were nontrivial true spatial 
covariances and no temporal correlation, while now an 
opposite situation will be considered — trivial true spatial 
covariances and a simple but fairly realistic model of tem- 
poral correlations, namely (cf. Sec. II B 2j) . the SVAR(l) 
with B = /SIjvj € (0,1), and idiosyncratic variances 
equal to each other. After a redefinition of parameters 
[into a dimcnsionlcss autocorrelation time, r = —1/ log/3, 
and variance a 2 = (c ld ') 2 / (1 — f3 2 )], the covariance func- 
tion (fT3|) is of the form (|102[) and reads 



C = 1 N , A(b-a) 



2 - 
a e 



(142) 



It is convenient to use the language of Fourier trans- 
forms (cf. Sec. Ill Dip , where (|76a|) . 



1 1 (a 1 
= — [A 2 -u 

A(u) Ai V u 



(143) 



where the two parameters of the model have been again 
changed to A\ = 2a 2 sinh(l/r) and A2 = 2 cosh(l/r), 
and the argument p to a variable on the centered unit 
circle [cf. (fT06|) ; v = u% 



, 1/2 1 r 

u ee e~ 2 ^, i.e., / dp(. ■■) = — 

-1/2 27T1 JC(0,1 



du-(...), 



u 



(144) 

in which it is most straightforward to apply the integra- 
tion method of residues. 
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FIG. 10: The critical rectangularity ratio (in the sector r < 1) 
above which the internal borderline arises for the TLCE for 
Toy Model 3. 



2. Master equations 



The discussion from Sec. Ill E 21 is relevant for this 
model — the master equations arc [ (|103ap - (|104p ] (with a 
complex unknown G and real nonnegative unknown h), 
which upon inserting (|143|) into them and using (| 144|) 
take on the form 



where for short, 



F^G, h)=rA\ 



Fi(G, h) =0, 
F 2 (G, h)=z, 



du- 



(145a) 
(145b) 



27ri7c(o,i) W(u) \G\ 2 + h' 



(146a) 



F 2 (G, h) = A 17 — 
Zm 



,2t-l 



du- 



(-u 2 + A 2 u- 1) 



C(o,i) 



W(u) 



(146b) 



where W(u) is a polynomial in u of order 2(i + 1) (i.e., 
4, 6, 8, etc., for t = 1, 2, 3, etc.), the leading coefficient 
{rA\G + Sti), and a property that its reciprocal coeffi- 
cients are mutually complex conjugate, 



W{u) 



rAiG - rAiA 2 Gu + rA x Gu 2 



2A- 



2 u 



t+2 



l A{[\G[ 

,t+3 



,t+l 



(147) 



2A 2 u ■ - + u~ 
rA\Gu 2t - rA 1 A 2 Gu 2t+1 



rAiGu 



21+2 



This polynomial has too high an order for it to be fac- 
torized analytically (thereby calculating the residues and 
consequently the two integrals), which is the main rea- 
son behind the difficulties in solving this set of integral 
equations — one has to resort to numerical methods. 



3. Borderline of the mean spectral domain 

The equation of the borderline is given by the above 
set with h = 0. Equation (|145a[) (a real equation for a 
complex unknown G = X + iY) gives thus a curve in the 
(X, Y)-plane, which is then mapped to the (x, y)-plane 
by Eq. (|145b[) . As evident from (|104p . the denominator 
(|147ll factorizes into two polynomials of order (t+1) each, 



W(u) 



(1 - A 2 u + u 2 +rA 1 Gu t+1 ) 



(148) 



(rAxG 



.t—i 



A 2 u l + u t+1 ) , 



whose roots (related by u — > 1 /u) should be found in or- 
der to evaluate the integrals [ (|146a)l . (|146b|) ] . 

Borderline for t = 1. For t = 1, one is faced with two 
quadratic equations, which therefore allows for consider- 
able analytical progress: The roots are easily derived, 
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FIG. 11: Monte Carlo eigenvalues versus theory for the TLCE (t = 1) for Toy Model 3. 



u 1± = A 2 /2±{Al/A-\-rA 1 G) 1 / 2 and u 2 ± = l/ul?, 
and so are the integrals in question. However, one has 
to carefully determine which two of these four roots lie 
inside the circle C(0, 1) and which two outside, for G sat- 
isfying (|145ap . A case-by-case analysis (performed in the 
most relevant situation of r < 1) reveals what follows: 

(i) The external borderline originates from Ui— and 
U2- inside C(0, 1). Its equation in the (X, y)-planc 
(|145ap reads 



ArA^X) U 2 ~ Ar AfY = 0, (149a) 



4 - A\ + ArAxX + 2r 2 A\ (X 2 + Y 2 



U 2 



-A\A 2 (X 2 + Y 2 ) U- 



l A 



Aj+ArAxX) X 2 - 



+ (-4 - A\ + ArAiX) Y : 



(149b) 



where U is an auxiliary real unknown (which can easily 
be solved out, which would however leave a very lengthy 
expression), and where the correct solution is such that 
the roots U2_ indeed lie inside C(0, 1). 
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FIG. 12: Monte Carlo eigenvalues versus theory for the TLCE (t = 5) for Toy Model 3. In (f), the inset presents the complete 
domain, while the main figure magnifies an area around zero. 



The map to the (x, y)-planc (|145b[l is 

1 



4 {X 2 + Y 2 ) ((1 + rA x X) 2 + r 2 A 2 Y 2 

((2/(rA 1 )+2X-rA 1 (X 2 + Y 2 )^jU 2 + 

h X (4 + Ai (-2(1 - r) + rA\) X- 
- 2r{l + r)A\X 2 



AA -2(l + r)+rAi-2r(l + r)A 1 X )Y Z , 



(150a) 



AY (X 2 + Y 2 ) ((1 + rA x X) 2 + r 2 A 2 Y 2 ^ 

- [x 2 (1 + tAiX) + (3 + tA x X) r 2 ) u 2 + 
X 2 (1 + rAxX) (-4 + A\ - ArA x X) + 
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FIG. 13: Monte Carlo eigenvalues versus theory for the TLCE (t = 15) for Toy Model 3. In (f), the lower inset presents the 
complete domain, while the upper inset and main figure magnify an area around zero. 



+ f - 4 + A\ + rAi (-12 + A 2 2 ) X+ 
+ 2r(l -3r)AlX 2 ^)Y 2 + 

+ 2r(l-r)A\Y 4 ^j. (150b) 

On the one hand, an algebraic form of these equations 
(in contrast to an integral form for t > 1) allows for quite 
easy plotting, as well as some analytical insight into the 
properties of the external borderline; on the other hand, 
any investigation of [ ([149ap - (|150b[l ] remains complicated. 



For instance, one may infer that the curve is symmetric 
with respect to the abscissa axis, in both coordinate sys- 
tems (Y — > —Y, which corresponds to y — > —y), and that 
it crosses this axis (Y = 0, which corresponds to y = 0) 
in two points (this number of crossings is a numerical 
observation and requires a proof); the crossings with the 
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X-axis are given by a quintic polynomial equation, 



(151) 



+ r 2 A\ {A% + 36r 2 - A 2 2 r 2 ) X*+ 

+ 16r 3 Al (8 - Af) X*+ 

+ 2r 2 A\ (4 - Af) (28 - A 2 2 ) X 2 + 

+ 12rA 1 (A-Alfx^ 

+ U-A 2 f = 



[one should select its real roots with the property that 
U2- lie inside C(0, 1)], which yield the crossings with 
the x-axis to be 



A- A 2 



2rA 1 X* + r(l - r)A 2 X 2 



2rA 1 (1 + rAM X 2 



[x+ may be checked to also obey a quintic polyno- 
mial equation, albeit lengthier than (|15ip . and thus not 
printed]. In Fig. [9l one sees plots of x* !m i n and x^ jmax 
versus r for r = 0.1,0.9, combined with the plots of the 
edges of the mean spectrum of the ETCE (cf. App. IC 3"|) . 
solving the sextic equation (|C7jl . [One can verify that — 
as may be suspected from Fig. [9}— in the limit St ^ 
(in which r ~ 1/St — > oo and r ~ 5t — > 0), the edges for 
both the ETCE and TLCE coincide; at the leading order, 
they both solve a certain quartic polynomial equation.] 
These graphs may be used as a genuine tool to deter- 
mine whether a given set of data has intrinsic temporal 
covariances, and of what r. 

(ii) The internal borderline arises when U2± are inside 
C(0, 1). It may be verified to exist only for 



1+C- 2 / 7 



< r < 1 



(153) 



(as mentioned, r > 1 has not been analyzed), which 
means that there occurs a "topological phase transition" 



(forming of a hole inside the domain) as one increases 
r past the critical value r c . This r c (r) is plotted in 
Fig- EH an d remark that it decreases from 1 to 1/2 when 
t grows from zero to infinity — this graph may also be a 
help in searches for true temporal covariances in multi- 
variate time series, yet only for r > 1/2. 

If (|153p holds, the equation of the internal curve in the 
(X, Y)-plane reads 



X 2 [ - 4 + A\ + 2(1 - 2r)A 1 X + r(l - r)A\X 2 

+ (a% + 2(1 - 2r)A 1 X + 2r(l - r)A 2 1 X 2S jY 2 + 
+ r(l - r)A\Y A = 0, 



(154) 



(152) for X_ < X < X + , with 



-l + 2r±y/l-r(l-r)Aj 

X± = -j . (155) 

7'(1 — r)A\ 

[The defining conditions of i*2± inside C(0, 1) are then 
automatically satisfied.] 

It is mapped to the (x, y)-plane by 

1 



(X 2 + Y 2 ) [p. + rAxXy + r 2 A\Y 2 
1 



2X + rA 1 (X 2 + Y 2 ) 

- X 2 (1 + rA^) (2 - A% + rA x X) + 

+ {Al + rA x (-3 + Al)X- 2r 2 A\X 2 ^Y 2 - 

- r 2 A 2 Y^j , (156a) 

-y 

(X 2 + Y 2 ) ((1 + rA x X) 2 + r 2 A\Y 2 
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2X + rA 1 (X 2 + Y 2 ) 



X(2 + rA x (1 + AT) X) + r A x (l + AT) Y 2 



(156b) 



Figure [TT] presents the numerically generated eigenval- 
ues of the TLCE and the numerically solved equations of 
the bo rderline(s), [ (fT49a|) - (fT50b)) ] and (for r > r c ) [ (|T54> 
()156b|) ] . for the parameters t = 1, and r — 0.1 (left col- 
umn) or r = 0.9 (right column), as well as r = 0.5 (top 
row), 2.5 (middle row) or 12.5 (bottom row), with very 
good agreement between the two, modulo the finite-size 
leaking out of the eigenvalues. The mean spectral do- 
main is "bullet-shaped" ; with growing r, it gets thinner 
in the y-direction (the more the smaller r) and longer in 
the positive part of the x-direction (the more the larger 
r). Moreover, Figs, [(e), (f)] show the phase with a hole 
inside the domain. 

Borderline for t > 1. In this case, one is left only 
with numerical methods to evaluate the borderline of the 
mean spectral domain, as the factors in (| 148|) have too 
high an order for an effective analytical treatment. The 
domain may have more complicated and very interest- 
ing shapes, possibly with several holes, as indicated by 
Figs. [12] and [T3] which correspond respectively to t = 5 
or t = 15, with other parameters as in Fig. [TT] 

All these graphs reveal that the spectrum of the TLCE 
is a sensitive tool for unveiling true temporal correlations 
present in the system, unlike the ETCE (cf. App. lC 3l and 
Fig. [24]), whose spectrum is "MP-shaped" and does not 
have enough conspicuous features to be an effective probe 
of autocorrelations. 

The numerical technique used to find these curves is 
very simple and could certainly be improved in order to 
obtain their more accurate approximations: (i) By trial 
and error, the extension of the borderline in the (X, Y)- 
plane has been estimated, (ii) The resulting rectangle in 
the (X, K)-plane has been divided into a grid of dimen- 
sions 1000 x 1000, and the value of F X (G, h = 0) (fTT5a|) 
at each point has been calculated by the residues, (iii) 
Only the values close to zero (with some assumed pre- 
cision) have been selected (|145a|) . (iv) These small val- 
ues have then been mapped to the (x, y)-values through 
[ (|145b|) . (|146b|) ]. yielding an approximation of the desired 
borderline. 

Remark (by inspecting Figs. [T2land [T3|) that the cen- 
tral region of the borderline [which may be checked to 
originate from its outer part in the (X, F)-plane] — which 
is also a region of the largest MSD (cf. Fig. fTT]) — is numer- 
ically reproduced much better than the peripheral region 
of the borderline [corresponding to its central part in the 
[X, y)-planc]. Hence, one may obtain a better approxi- 
mation to the borderline by using an (X, Y)-grid denser 
in the central part. 



4- Mean spectral density 

Calculating the MSD is even more involved than the 
borderline as one needs to solve the master equations 
[ (I145all - (|145b|) ] for both unknowns, G and h — the nu- 
merical algorithm is outlined in App. [D] and Fig. [TT] vi- 
sualizes the resulting MSD in 3D. 



E. Toward more realistic toy models 

In the above examples, the effects of either the spa- 
tial or temporal intrinsic covariances on the MSD of the 
TLCE have been separately analyzed. However, in a de- 
scription of complex systems such as financial markets 
one will certainly encounter spatial and temporal covari- 
ances both present and mixed (cf. Add. I A 2"j) . such as in 
the VAR models (cf. Sec. II B 2[) . Such situations pose 
a much greater challenge for an analytical approach to 
the MSD, and another publication should be devoted to 
them, but in this Section, some hints will be given as to 
the directions of generalizations, as well as a few prelimi- 
nary results. Moreover, the ETCE will also be discussed 
since its MSD for the below models have as yet been 
unknown. 



1. Toy Model J^a: Different variances and identical 
temporal exponential decay 

Definition. A first step toward more meaningful mod- 
els could be to combine the above Toy Models 2a and 3, 
i.e., suppose that the assets are divided into K sectors of 
distinct variances a 2 (hence, the different assets are still 
independent), all of them exhibiting exponentially decay- 
ing temporal autocorrelations, with the same autocorre- 
lation time r. This is the SVAR(l) model (cf. Sec. [TB2]) 
with B proportional to the identity matrix and a finite 
number of distinct idiosyncratic variances. The covari- 
ance function (|12p is factorizcd and translationally sym- 
metric in time (1831) with 



C = diag( a\ 



>K 



A(b-a) 



N\ times N2 times 



Nk times 



(157) 



ETCE. The master equation for M = M c (z) is ([4"8l 
with (fT34|) and (|C5|) . that is 



K 



Pk 



k=l ro%M(x+s) 
2 1.9 



2 M 2 



r 



X 



1. 



(158a) 
(158b) 



where s is an auxiliary complex unknown, and recall 
Pk = Nk/N, x — coth(l/r). By solving out s, one ob- 
tains a (lengthy) polynomial equation for M of order 
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FIG. 15: Monte Carlo eigenvalues versus theory for the ETCE for Toy Model 4a. The graphs in the main figures are hardly 
distinguishable, hence the insets show them separately. 
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FIG. 16: Monte Carlo eigenvalues for the TLCE for Toy Model 4a. 



2(K + 1). which is easily appropriated numerically. [Re- 
mark that switching on exponentially-decaying autocor- 
relations doubles the polynomial order of the master 
equation for the ETCE; compare the current situation 
with Toy Model 2a JC4|, and Toy Model 3 |C6} with 
Toy Model 1 rtCil .] 

In Fig. Ql)l the MSD from Monte Carlo simulations is 



compared to the numerical solutions of [ (|158ap . (|158b|) ] . 
for the values of the parameters K = 5, erf = 1, a\ = 2, 
0-3=4, cr| = 6, erf = 8, equal relative multiplicities 0.2, 
and three values of r = 0.02,0.1,0.5, all for cither r = 
2.5 (a) or t = 12.5 (b), and all with perfect agreement. 
One observes that the ETCE could still be an efficient 
tool in distinguishing the variance sectors (the distinct 
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"lumps"), the better the smaller r, although the effect 
of the autocorrelations is to smear these lumps out, the 
more severely the greater r, making this analysis more 
difficult. 

TLCE. The master equations arc [ (|95a[) - (|95c[) ]. with 
POaj) . (pObf . flSSJ}, (HSU, (|94b)) . ([52])]. and a numerical 
method of solution remains to be developed. 

Figure [16] only presents the numerical eigenvalues of 
the TLCE, for t = 1, K = 2, of = 1, of = 5, equal rel- 
ative multiplicities 0.5, as well as r = 0.1 (left column) 
or r = 0.9 (right column), and r = 2.5 (top row) or 
t = 12.5 (bottom row). The domain is "bullet-shaped," 
with a possible hole inside, similarly as for Toy Model 3. 
However, for sufficiently small r and r, one recognizes its 
parts corresponding to the different variances; this shape 
might be exploited in determining the parameters of the 
model from a given set of data, once the theoretical mas- 
ter equations become numerically solvable. Also, just like 
for the ETCE, increasing r or t makes these parts less 
distinguishable. 



2. Toy Model 4b: Different variances and different 
temporal exponential decays 

Definition. A reasonable way to further generalize the 
above model would be to assume that each variance sec- 
tor has also a different autocorrelation time, 



C(c) = diag( ^^'^f^_^> 

JVi times N2 times 



' K K 

Nk times 



(159) 



This is the SVAR(l) model with B diagonal with a fi- 
nite number of distinct beta coefficients corresponding 
to distinct idiosyncratic variances. The covariance func- 
tion (|12p is translationally symmetric in time but not 
factorized into a spatial and temporal parts, being the 
first instance of such mixing in this paper (it is Case 3; 
cf. Sec. ED). 



The Fourier transform (|76a[) of (|159[) in the variable u 
flUUl reads [cf. dT43ll ]. 



C(u) = diag 



A 2 ,i 



A 



2,K 



N± times 



N K ti: 



(160) 

with Ai k = 2of sinh(l/Tfc) and A^^k = 2cosh(l/rfe). 

ETCE. There is a common pattern in solving the mas- 
ter equations (for both the ETCE and TLCE) in the sit- 
uation of mixed spatial and temporal covariances, which 
will be sketched on the current example: Firstly, the 
matrix structure of C(u) determines the matrix struc- 
ture of E™, as visible from (|77a|> . i.e., S** is diago- 
nal, with some elements with multiplicities Nk, for 
k = l,2,...,K; these constitute the unknowns. Sec- 
ondly, (|78ap implies that the same is true of G NN , with 
G k = l/(z - E fe ). Thirdly, (f77bj) becomes 



(«) 



K 

1=1 



l 



A 2 ,i -u - 7; z-Ei 



(161) 



Fourthly, one plugs this into (|78bp . and fifthly, further 
into (I77a|) . which gives the final K equations for the K 
unknowns, 



1 

2?ri 



C(0,1) 



du 

u A 



Ai„ 



1 



2,k 



u — - 



E tt (m) 



(162) 



Solving this system requires computing the residues of 
the integrands inside the centered unit circle, which is 
equivalent to finding the zeros in u of their denominators, 
and these are all polynomials of order 2K; therefore, it 
poses a challenge to construct even a numerical algorithm 
able to calculate the values from (|162|) . 

In Fig. [T71 there are examples of Monte Carlo his- 
tograms of the MSD of the ETCE, for K = 2, of = 1, 
of = 5, equal relative multiplicities 0.5, three values of 
r = 0.02,0.1,0.5, and respectively T\ = 0.5, t-2 = 2.5 
(a) or n = 2.5, r 2 = 0.5 (b). As with Toy Model 4a 
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FIG. 18: Monte Carlo eigenvalues for the TLCE for Toy Model 4b. 



(cf. Fig. IT5j) . the variance sectors are visible for small 
enough r, but larger autocorrelation time broadens 
the peak corresponding to <7&. 

TLCE. The master equations are [([751 ]) - ([52 ]) ]. and a 
similar method as for the ETCE should be used to re- 
cast them in a form appropriate for a yet to be devised 
numerical procedure. 

In Fig. [TU the Monte Carlo eigenvalues of the TLCE 
are shown, for t = 1, K = 2, of = 1, erf = 5, equal rela- 
tive multiplicities 0.5, and r = 0.1, as well as six combi- 
nations of the values 0.5, 2.5, 12.5 for n, T2- One again 
recognizes that increasing any Tk causes the sub-domain 
corresponding to a\ to expand, with rate depending on 



the value of the variance (it seems it has a stronger im- 
pact on smaller variances) ; for sufficiently small autocor- 
relation times (especially for smaller variances), and of 
course for small enough r, the sub-domains may sepa- 
rate (topological phase transition), which has not been 
observed for Toy Model 4a. 



3. Toy Model 4c: A simple vector autoregression model 
with a market mode 

Definition. Continuing the above progression of 
SVAR(l) models, a next reasonable step could be to in- 
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FIG. 19: Monte Carlo eigenvalues versus theory for the ETCE for Toy Model 4c. The insets magnify a small peak originating 
from the "large" eigenvalue As(c), while the peak around zero corresponds from the two "small" eigenvalues Ai,a(c). 



elude genuine correlations between various assets. A sim- 
plest such example is to consider the market mode which 
depends on itself in the previous moment of time, while 
any other asset depends both on itself and the market in 
the previous moment, as outlined at the end of Sec. lIB2l 
This is the most advanced model in this article, with the 
time-dependent covariance matrix [ (|14a|) - (|14e|) ] (belong- 
ing to Case 3). 

The Fourier transform (|76a[) of the covariance function 
(|TT|) of an arbitrary VAR model reads 



C(p) 



1 



K(p) 



1 



1 - K(p)t 



(163) 



which in the considered simplified version turns into 
(note, K(u) = uB), 

Cn(u) = ci{u) = 

u (- 7 m 2 + (1 + 7 2 + (N - l)f3 2 ) u - 7) 



Cu{u) 
C u (u) 
C a (u) 



(u — a)(ua — l)(u — 7)^7 — 1) 
— a 

cs(u) = 



(u-7)(u 7 - 1)' 
-f3u 



c 4 (it) 



(u — a)(u — j)(wy — 1) ' 

^ 

(ua — l){u — 7) (1*7 — 1) 



(164a) 
(164b) 

(164c) 



(164d) 



for i > 1, with the remaining matrix elements zero (i.e., 
there are four distinct matrix entries on the diagonal, first 
row and first column). This constitutes an input data 
for the master equations [ (jT?aj) - (|78bj) ] (for the ETCE) or 
[ ([79ajl -p])] (for the TLCE). 

ETCE. Following the same logic as in Sec. IIIIE21 
one starts from noting that according to (|77ap the ma- 
trix structure of S™ mimics that of C(u), i.e., it 
is determined by four complex parameters, = £1, 



^NN 



£2, £1 



NN 



£ 3 , £ 



NN 



£4, for all i > 1, which 



are the unknowns. One substitutes this form into (I78al) . 



where the matrix inversion is explicitly doable, and fur- 
ther into (|77bj) . 

£ TT (u) = 



S3S4 



( z - s i)at=t 



z-S 2 

c 2 (u)£ 3 £4 N - 2 
(z-£ 2 ) 2 N 
c 2 (u)(z - £1) + c 4 (m)£ 3 + c 3 (m)£ 4 1 



z — £9 



N 



Ci (u) 



1 



N(N — 1) 



(165) 



[Although N is always infinite in this paper (O , here one 
should treat it as finite, governing the number of assets 
relative to the market size.] Plugging this result into 
(|78b[) . and then into (|77a[) . one arrives at four integral 
equations for the four unknowns, 



1 

27ri 



du— c s (u)- 
C(o,i) u 



T, TT (u)' 



(166) 



for s = 1,2,3,4. One finds that the four integrands are 
rational functions of u, with the same denominator being 
a quartic polynomial in u. For that reason, analytical 
solution is impossible, and an effective numerical method 
needs to be devised. 

Figure [19] presents the Monte Carlo histograms of the 
eigenvalues of the ETCE, for r = 0.1, a = 0.9 and respec- 
tively /3 = 0.1, 7 = 0.8 (a) or /3 = 0.8, 7 = 0.1 (b). The 
eigenvalues Ai i 2, 3 (c = 0) [ (fl~5|) , (fl"tj|) ] of the true covari- 
ance matrix acquire then the approximate values 2.7778, 
2.1005, 94.8502 (a) or 1.0101, 1.0082, 408.735 (b), which 
are seen to be broadened by the measurement noise in the 
form of a large peak corresponding to the two smallest 
eigenvalues (they arc close to each other and thus in- 
distinguishable from this spectrum) and a small distant 
peak originating from the one large eigenvalue. 
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FIG. 20: Monte Carlo eigenvalues for the TLCE for Toy Model 4c. The main figures present the complete spectrum, while 
the insets zoom into the central region, corresponding to Ai,2(c). There are 100 (= number of iterations) points in the distant 
region, corresponding to As(c). 



TLCE. The master equations in a form suitable for nu- 
merical evaluation are derived analogously to the ETCE 
case above; they will not be explicitly written down, only 
some examples of the Monte Carlo spectra of the TLCE 
are depicted in Fig. [201 for the same values of r and the 
beta coefficients as above, as well as respectively t = 1 
(top row), t = 5 (middle row) or t = 15 (bottom row). 
Again, there is a bullet-shaped central region correspond- 
ing to the two smallest eigenvalues and a distant region 
from the noise-broadened large one. In other words, one 



not only learns from these plots about the presence of 
the market mode, but by investigating the shapes of the 
domain for various t one may infer about the autocorre- 
lations in the system. 
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IV. CONCLUSIONS 
A. Summary 

The focus of this article is on an important problem of 
assessing from historical data correlations between var- 
ious objects of a complex system (such 8jS 0, financial 
market) at identical or distinct time moments. In the 
first case, a useful tool is the spectral analysis of the 
equal-time covariancc estimator, which is a Hcrmitian 
random matrix and which has been extensively studied 
in the literature. It has been shown for various mod- 
els to reproduce the true equal-time correlations, albeit 
with a measurement noise, which nevertheless diminishes 
with decreasing rectangularity ratio r. Unfortunately but 
quite understandably, its spectrum practically fails to be 
sensitive enough to true time-lagged correlations, a prop- 
erty perhaps even more important, reflecting rich, com- 
plicated and relatively little known temporal dynamics 
of the markets (e.g., heterosccdasticity, leverage effects, 
Epps effect etc.). 

It is the time-lagged covariance estimator which con- 
stitutes a sensitive probe of these temporal correlations. 
It is a non-Hermitian random matrix, and as such is 
much more involved from the analytical point of view 
than Hcrmitian matrices, which explains why it has been 
largely neglected. The only analytical result to date has 
been [4lL |42| , whose main result is shown here to be in- 
correct while the right one is given. 

This paper is a comprehensive treatment of both these 
estimators. The method of planar diagrammatic expan- 
sion, adjusted to the non-Hermitian world (explained at 
length in App. [Bj) , is applied in the most general Gaus- 
sian case to derive the master equations (cf. Sec. |TT] and 
Tab. H| from which the mean spectral density can be cal- 
culated. 

These general results are then used for a detailed spec- 
tral analysis (cf. Sec. IIIII and Tab. Ill|) — i.e., finding the 
density, comparing it to Monte Carlo simulations and 
discussing it with respect to the effectiveness in unveil- 
ing true covariances — of the time-lagged (and sometimes 
also equal-time) covariance estimator for four toy mod- 
els of the market, each designed to approximate its cer- 
tain aspect. The simplest Toy Model 1 has zero true 
correlations and may be treated as the null hypothesis. 
Toy Model 2 assumes zero temporal correlations and in- 
vestigates nontrivial variances, either their finite number 
(2a) or a power-law distribution (2b), which both model 
the presence of industrial sectors. They are recognized 
to be much more accurately reflected in the spectrum 
of the equal-time estimator. Toy Model 3, on the con- 
trary, has zero correlations between distinct assets but 
nontrivial autocorrelations; it is a simplest example of the 
class of vector autorcgrcssion models, with exponentially- 
decaying temporal correlations identical for all the as- 
sets. The spectral analysis is quite difficult both analyt- 
ically and numerically but leads to very intricate shapes 
of mean spectral domains which may serve as a sensitive 



probe of true temporal correlations in the system, unlike 
the equal-time estimator. 

Since these models are both already complicated and 
yet quite crude approximations, three more realistic mod- 
els are further investigated — due to analytical difficul- 
ties, both the equal-time and time-lagged estimators are 
considered but it is only sketched how to write down 
their master equation in a form suitable for numerical 
evaluation and Monte Carlo simulations are presented. 
Toy Model 4a has a finite number of variances and 
exponentially-decaying autocorrelations identical for all 
the assets; Toy Model 4b differs from it by assuming dis- 
tinct autocorrelation times for all the sectors; Toy Model 
4c is the most advanced example with a nondiagonal 
time-dependent covariance matrix originating from the 
presence of the market mode such that the evolution of 
each asset depends on the value of this asset and of the 
market in a previous moment of time. 

Due to analytical tractability, all the above calcula- 
tions are based on the Gaussian probability distribution 
of the assets, which is an assumption universally made 
both within the scientific approach and financial prac- 
tice, yet certainly quite far from the truth (cf. App. [A"l 
however, this can be a basis for non-Gaussian generaliza- 
tions through random volatility models). Nevertheless, 
four extensions of Toy Model 1 beyond this realm are also 
accomplished: (i) Gaussian assets but with the EWMA 
estimator; (ii) Student t distribution of the assets, ob- 
tained either from a common random volatility factor or, 
more realistically, IID temporal volatilities; (iii) free Levy 
distribution, which is a more tractable approximation to 
the more realistic Wigner-Levy distribution of the assets. 
It is the first step in describing how fat tails impact the 
spectrum of the time-lagged estimator; applications to 
more complicated toy models constitute a task for an- 
other project. 

An addition to the above results is a formulation and 
numerical tests of two mathematical hypotheses valid 
for non-Hcrmitian random matrix models with rotational 
symmetry in their mean spectrum: (i) the TV-transform 
conjecture, which relates the eigenvalues and singular val- 
ues, thereby being a bridge between Hermitian and non- 
Hermitian worlds; (ii) the form-factor using the comple- 
mentary error function which describes the universal de- 
cline of the density close to the borderline of its domain. 



B. Open problems 

Several paths may be taken from where this paper 
ends. One major research program would be to repeat 
the analysis of the mean spectral density for more realis- 
tic (say, from the financial point of view) toy models; this 
may include: (i) Development of numerical methods able 
to handle the master equations in the Gaussian case but 
for e.g. the vector autorcgrcssion (VAR) models with 
more involved structure of cross-covariances; the "vec- 
tor moving average" (VMA) models; or a combination 
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of the two called the "VARMA" [U . (ii) Generalization 
to weighted estimators, such as the EWMA or with log- 
arithmic weights, (iii) Generalization to non-Gaussian 
monovariate probability distributions of the assets (e.g., 
Student or Levy) for models with nontrivial covariancc 
structure, (iv) Attempts of incorporating more involved 
temporal effects, such as the heteroscedasticity, mono- 
variate or index leverage effects, Epps effect etc. 

An important part yet to be accomplished is a com- 
parison of the results stemming from some selected toy 
models with real-world financial data, and learning how 
to validate them and assess their parameters by fitting. 
This method may give a profound insight into the struc- 
ture of intrinsic correlations between financial objects. 

Besides, several mathematical problems are worth in- 
vestigating with a view for a broader understanding of 
the time-lagged estimator; one could compute: (i) finite- 
size effects; (ii) statistical properties of the eigenvectors; 
(iii) the statistics of the "largest" (most distant) eigen- 
value, such as in Toy Model 4c; (iv) a more detailed 
analysis of the topological phase transitions in the mean 
spectral domains. 
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Appendix A: Non-Gaussian probability distributions 
for financial instruments 

This appendix is a short outline of the main features 
of the probability distributions encountered in finances; 
they are non-Gaussian and reflect an involved structure 
of correlations present on the market. On the contrary, 
in the majority of this paper the Gaussian distribu- 
tion is assumed — with notable exceptions of Sees. MI B 71 
(Student) and IIIIB 81 (free Levy) — with relatively sim- 
ple structure of the true correlations. This is because 
the calculations in this approximation are already very 
complicated, yet much less than with any of the addi- 
tional features taken into account; also, Gaussian results 
can serve as a basis to construct more realistic random 
volatility models (cf. Sec. IA 2[) . This appendix explains 
what sorts of errors are committed by using the Gaussian 
distribution, and the necessary means to amend them. 



For a more thorough introduction to the subject, cf. es- 
pecially the textbook [52| , as well as the reviews 11. 281. 

1. Student t and free Levy distributions 

a. Weak stability for covariances 

In the statistical approach, one typically makes an 
assumption of "weak stability," i.e., that the mecha- 
nisms governing the market (thus, the probability dis- 
tributions) do not significantly vary with time over the 
considered periods, so that historical data could be used 
to assess future behavior (cf. Sec. II A 2|) . This is rather 
untrue during the crisis times, and also in attempts to 
predict future returns (which have highly non-stationary 
distributions), but is typically justified for covariances 
(thus, risks), which are the subject of this article. There- 
fore, the mean trend ("drift") will be subtracted in all the 
considerations and fluctuations will be the only focus, 

(Ria) = 0. (Al) 

A further comment on these disregarded mean values 
of the returns: First of all, for short time periods (less 
than several weeks) , the mean return is negligible as com- 
pared to its "volatility" (a financial term for the stan- 
dard deviation; e.g., for stocks over one business day, 
a fraction of a percent versus several percent). More- 
over, the returns should have no temporal autocorrela- 
tions, since their nonzero value would enable a profitable 
trading strategy ("arbitrage") until these autocorrela- 
tions would disappear, according to the "efficient market 
hypothesis." On the other hand, there exist riskless as- 
sets (bonds) , hence stocks should earn on average at least 
the riskless rate of return; furthermore, there are factors 
impacting the mean returns which are not random and 
quite predictable, such as short-term interest rates. For 
these reasons, the autocorrelations of returns are actu- 
ally nonzero, yet weak and over small time lags (less than 
about 30 minutes for liquid assets). Note that this does 
not contradict the efficient market hypothesis — as it may 
seem, profits could be obtained using the high-frequency 
(HF) trading — because of transaction costs. They also 
have an important implication, being one reason for the 
"Epps effect" (cf. the end of Sec. IA 2 b|) . Some risk as- 
sessment methodologies (such as RiskMctrics 2006 (53|) 
attempt to incorporate these effects; however, they will 
henceforth be neglected. 



b. Fat tails 

In order to make the discussion from this paper truly 
practically applicable, one will need to face the fact that 
the returns are strongly non-Gaussian. It is an observed 
fact that their unconditional probability distribution over 
all the time scales between several minutes and several 
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days has "fat (heavy/Paretian/power-law) tails," 

1 



unconditional 



\R.. 



, R„ -> ±00, (A2) 



where the experimental value of the exponent is fx ~ 
3 -T- 5 [U, [55| . This reflects considerable chances for ex- 
treme events, and also signifies the presence of multi-scale 
phenomena (i.e., when both very small and very large val- 
ues appear). For any such distribution, all the moments 
(i2?.) for n > fx diverge. 



c. Student t distribution 

The financial data can be well fitted with the Student 
t-distribution (cf. e.g. HUSUl), 

^unconditional = ~~J= _ / „ \ T+JT 1 (A-3) 

r (f) (6 2 + R 2 „) — 

where T(-) here is the gamma function, and the variance 
reads a 2 = 9 2 /(fi - 2). 



d. Free Levy distribution 

Wigner-Levy random matrices. Another candidate 
for a fat-tailed probability distribution of the returns 
is the Levy stable law. Recall that it is described by 
four parameters, the stability index a £ (0,2], skew- 
ness (3 £ [—1,1], range 7 > and location 5 £ R 
(which will always besides this Section be set to zero). 
Its characteristic function (related to the PDF through 
L a ,p,-y,s(x) = ^ J**™ dkc~ lkx L a .fj tl . s (k)) reads 

log L a ,p,7,s(k) = 

( -7|fc| Q (1 -sign(fc)i/3tan^) + iSk, for a ^ 1, 
~ \ -j\k\ (1 + sign(fc)i/3f log [fe|) + iSk, for a = 1. 

(A4) 

It has: (i) fat tails, L a ^ n ^{x) ~ ^4±/|£| 1+Q for x — > ±00; 
(ii) all the positive moments divergent for a < 1, and 
only the mean finite (= <5) for a > 1; (hi) the relative 
"weight" of the right and left tail quantified by skewness, 
(A+ - A-)/{A+ +A-) = P; (iv) the scale 7 1 /". Besides, 
it possesses a host of mathematically interesting features, 
in particular: (i) It is a fixed point of convolution mod- 
ulo an affine transformation (dilation and translation), 
i.e., the PDF of a sum of a large number T of IID ran- 
dom variables is identical to the PDF of the single vari- 
able modulo an affine transformation if and only if the 
variables are Gaussian or Levy. The dilation is "super- 
diffusive," i.e., this sum should be rescaled by T 1-1 /". (ii) 
It is a-stable, i.e., a sum of a-Levy (not necessarily inde- 
pendent) remains a-Lcvy. (iii) It constitutes the basin of 
attraction for summing a large number T of random vari- 
ables (not necessarily independent or identical), rescaled 



super-diffusively, which have fat tails with the exponent 
in (0,2] [a generalized central limit theorem (CLT)]. 

In RMT one wants to mimic the above scalar construc- 
tion while preserv ing its foundational properties. The 
most obvious way [56| (cf. also [Ht], HH) would be to fill a 
real symmetric T xT matrix with IID (for a < b) entries 
h a b sampled from some Levy stable distribution, rescaled 
super-diffusively to ensure a finite thermodynamic limit, 



ttWL 
H ab 



h a b 
yl/a ' 



(A5) 



which is dubbed a "Wigner-Levy" (WL) or "Bouchaud- 
Cizcau" matrix. This may then be extended to complex 
or rectangular matrices. 

Free Levy random matrices. There exists however an- 
other matrix generalization of the Levy stable laws in 
which one asks about matrix ensembles whose MSD is 
preserved modulo an affine transformation upon sum- 
ming IID random matrices. To be precise, the notion 
of statistical independence has been generalized to the 
noncommuting realm, and called freeness, in Voiculescu- 
Speicher free probability theory [2^, |3(|. Moreover, in 
analogy to the prescription in classical probability the- 
ory for summing IID random variables (one should sum 
the logarithms of the characteristic functions), in order 
to sum free Hermitian random matrices Hi and H2, one 
should invert functionally their holomorphic Green func- 
tions, 

Gh(BhW)=Bh(GhW) = z, (A6) 

which is called the "Blue function" and which in turn 
obeys the "addition law," 

1 



-Bh!+h 2 (z) = B Hl (z) + Bh 2 (z) - 



(A7) 



[Note a close relationship with the multiplication law 
(l49|) . outlined in Sec. Ill C Tl which uses the TV-transform 
(|47]) being the functional inverse of the holomorphic M- 
transform. Cf. also the end of App. IB 2 al for its non- 
Hermitian generalization.] 

This new stability condition can then be solved with 
aid of (|A7|> leading to the "free Levy stable laws" [H| , 



f S + bz 01 " 1 



for a^l, 



a-i 7 (l + 0)-2£*log(7*) + i for 



1, 



(A8) 



where the branch structure is chosen such that the upper 
complex half-plane is mapped to itself, and where for 
short, 

r 7 e^W», fora<l, 

^{^(H^), fora>1 . (A9) 

The free Levy laws are convenient to work with due to 
their free probability framework, and actually do not sub- 
stantially differ from the WL distributions [58|. There- 
fore, in a series of works (cf. e.g. j3ll l60l - [63 |). they have 
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been chosen to model the behavior of financial assets. 
E.g. in [U, the MSD of the ETCE observed for TV = 406 
stocks from the S&P 500 market over T = 1309 days (i.e., 
r ~ 0.31) has been well fitted with the free Levy distri- 
bution with a = 3/2 and j3 = 0. 



2. Random volatility models 

a. Monovariate random volatility models 

However, the above no n- Gaussian probability laws are 
only a part of the picture, since it is well known that 
the returns at different time moments are not IID, but 
have different distributions and are correlated with each 
other. Otherwise, their sum would too quickly converge 
toward a Gaussian fixed point (according to the CLT for 
IID variables) , contrary to the observation (vide: market 
crashes). A convenient way to express this is in terms of 
"random (stochastic) volatility models" [HH (a.k.a. "de- 
formations" [49l. I65l468l |). where the return at time a, 



R a 



(A10) 



(the index i is removed in this Section), where the stan- 
dardized returns ( "residuals" ) e a arc assumed IID Gaus- 
sian N(0,1), while the volatilities a a are also taken to 
be random, where the JPDF V({e a }, {er a }) should repro- 
duce the following "stylized facts" : 

(i) It should yield the proper marginal distribution of 
the returns. For instance, if all the a a were independent 
from each other and from all the e a (which is not true, but 
may be considered a simple toy model), then choosing 



1 



6"' 



2/V2-ir 



(All) 



(i.e., the inverse variance has the gamma distribu- 
tion, which is a two-parameter generalization of the 
chi-squarcd distribution) would lead to the Student t- 
distribution ([A"3]) , 

(ii) The volatility depends on time (in a random way) 
and displays time-lagged correlations, which are quite 
small (a few percent; this value depends on the volatil- 
ity proxy used — e.g., the "HE proxy," being the daily 
average of HF volatilities — and the time lag), but per- 
sisting over long periods of time ("long memory") — i.e., 
that high volatilities persist in crisis times and low ones 
in calm periods — which may be described by a power-law 
(slow) decay, 



= a(\b- 



(A12) 



where the averaging uses the volatility JPDF, and the 
exponent v ~ 0.2 0.4. This phenomenon is called 
"heteroscedasticity," "volatility clustering" or "intcrmit- 
tence," in analogy with the motion of liquids, where also 



periods of laminar and turbulent flow persist consecu- 
tively (cf e.g. miMIzii)- 

Remark 1: In particular, it implies that estimation 
based on historical time series is a valid means to assess 
future risks. 

Remark 2: Even if the residuals were uncorrelated with 
the volatilities (which is not true; see below), then the 
returns would also be uncorrelated, (R a Rb) = 5 a b{cr a (Jb), 
but they would not be independent, hence there would 
be nontrivial higher-order correlations, such as (R^Rf) — 
(Rl)(Rl) = g(\b-a\), for a ^ b (from the Wick theo- 
rem). 

Remark 3: With the heteroscedasticity (|A12|) present, 
a sum of a large number T of returns (aggregation with 
time) still tends to a Gaussian fixed point (CLT), but 
much slower than for independent variables. One may 
check this by computing the kurtosis kt (the normalized 
fourth moment; it measures a distance of a given distri- 
bution from the Gaussian) of this sum, which for large T 
behaves as 




for v > 1, 
for v < 1, 



(A13) 



which means that the convergence is slow (several weeks) 
for long-memory volatility autocorrelations ( "anomalous 
kurtosis"). 

Remark 4 : A common approach to modeling the 
volatility dynamics is by a stochastic process. Even 
though a proper way would be to search for models of 
the mechanisms underlying the volatility process, one of- 
ten selects a process according to the criterium of com- 
putability, which in practice restricts their class to the 
so-called (integrated) "auto-regressive conditional het- 
eroscedasticity" (ARCH) processes [7114731 . in which the 
variance depends linearly on the past squared returns, 



T 

£ 

5=1 



WbR a -b 



(A14) 



where T here is a lower cutoff, and Wb are some posi- 
tive weights, satisfying the sum rule, J2b=i w b = !• An 
example broadly used by the financial industry (intro- 
duced in RiskMetrics 1994 ||, Q) was the EWMA flU 
(cf. e.g. [ii, (zH), Wb oc exp(— b/r), where one parame- 
ter is used for all the assets, rdt = 16.2 business days. 
In RiskMetrics 2006, it was argued that a better fit is 
given by a logarithmic decay, Wb oc 1 — log bj log r, with 
rSt ~ 3 -j- 6 years. 

(iii) The "leverage effect" [76T - [78j . i.e., the empirical 
fact that negative past returns typically increase future 
volatilities, which means that e a and a a are correlated in 
such a way that 



< 0, for b > 



(A15) 
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b. Multivariate random volatility models 

The above experimental constraints on the statistics 
of a single return are both very nontrivial to implement 
in such a way that the resulting model would be mathe- 
matically tractable, and also leave plenty of freedom for 
model specification. To make things even tougher, this 
is still an incomplete description, as there are many cor- 
related assets on the market, 



Ria 0~ia£ia 



(A16) 



("multivariate random volatility model"), and one has 
to select a JPDF which will additionally reproduce their 
observed dependence structure. 

Multivariate Gaussian returns. This task is straight- 
forward only if the returns are Gaussian (which, as ex- 
plained above, is quite far from truth), since the most 
general (real) multivariate Gaussian distribution (of zero 
mean) is fully described by the two-point covariance func- 
tion, 



^ij,ab — (RiaR 



(A17) 



Pseudo-elliptical and elliptical models; multivariate 
Student returns. This suggests a way to introduce a de- 
pendence structure between the assets using the form 
(IA16|) (the index a is removed in the following discussion, 
as it pertains to the cross-correlations only) — by sup- 
posing that the Gaussian residuals are cross-correlated, 
(eiCj) — Cij, and multiplying them with some random 
volatilities (Xj, perhaps also cross-correlated, but indepen- 
dent of the residuals — which leads to the class of "pseudo- 
elliptical models." 

A commonly used example ("elliptical models") is to 
take a single volatility for all the assets, Ri = ere;. For 
instance, if the volatility is distributed according to ([Alip 
(with 2 = 2/i), this yields the "multivariate Student t- 
distribution" of the (simultaneous) returns, 



V{{Ri}) = 



v^pDctc r(f) 
1 



(A18) 



In this case, the matrix C is sufficient to describe all 
the correlation functions: the two-point one, (RiRj) = 
j^Cij, as well as the higher ones th roug h the Wick the- 
orem. However, it has been argued [79| that the JPDF 
of returns is not elliptical, and several volatility factors 
need to be exploited. 

Copulas. Let it also be added that for more gen- 
eral non-Gaussian distributions, the two-point covariance 
function is typically insufficient (because of higher-order 
correlations), invalid (because it diverges), or uninter- 
esting (e.g., because one is more concerned about cor- 
relations of extreme events). There is a method — based 



on the Sklar theorem [80J — which factorizes the problem 
of defining multivariate non-Gaussian distributions into 
a specification of marginal distributions [such as (|A3|) ] 
and a description of dependence: (i) One transforms 
the returns Ri into Ui = V< t i(Ri), where T J < ^(-) is the 
marginal cumulative distribution function (CDF) of asset 
i; these new random variables are uniformly distributed 
on [0,1]. (ii) Their joint CDF, 



C({w i })=Prob({[/ i <'U i }): 



(A19) 



is named the "copula," and is responsible for the depen- 
dence of the returns. This framework includes a number 
of known measures of dependence used in the presence of 
fat tails, such as the Kendall tau, Spearman rho or tail 
dependence. A copula has to be chosen so as to repro- 
duce the empirical two-point covariance function (RiRj), 
which again is a greatly underconstrained task ("cop- 
ula specification problem" ) , leading to the choices based 
rather on mathematical computability, and not founded 
upon a solid financial reality (cf. e.g. (8ll. l82j). 

Temporal cross-correlations. From the above dis- 
cussion, it is quite obvious that disregarding temporal 
correlations — not only the autocorrelations of any sin- 
gle asset but between distinct assets at different time 
moments — is a serious simplification. There are at least 
two important stylized facts which require incorporation 
of the temporal cross-correlations: 

(i) The "index (correlation) leverage effect" (83l-[85|— 
not only do the volatilities increase upon large negative 
past returns [cf. point (iii) in Sec. I A 2 a| . but so do the 
cross-correlations ("cross-correlations jump to one in cri- 
sis periods"). This is reflected by the fact that the lever- 
age correlation (|A15[I for market indices (i.e., the aver- 
ages of the returns of all the stocks) is stronger than for 
separate stocks — the additional contribution originating 
precisely from the cross-correlations. 

(ii) The "Epps effect" H(| — the equal-time cross- 
covariance (Ri a Rja) grows significantly as one increases 
the observation time lag St between about 5 and 30 min- 
utes, then more slowly, finally to saturate after a few 
days (cf. also (87|)- It can be explained [H, \8^, by 
relating (Ri a Rja) to the time-lagged cross-covarianccs 
(RiaRj.a+t) on shorter time scales 8t, and proposing rea- 
sonable toy models of the latter, which reproduces the 
"Epps curve," in agreement with the experiment. It is an 
important feature of the structure of cross-covariances, 
and also from the point of view of their historical mea- 
surement (cf. Sec. II A2l . 



Appendix B: Gaussian planar diagrammatic 
expansion 

This appendix contains (i) an introduction of the basic 
notions of random matrix theory in whose language the 
results of this paper are expressed (Hcrmitian case in 
App. IB 1 a| non- Hcrmitian in App. IB 2 a|) : (ii) a sketch 
of the main technique used to obtain these results, i.e., 
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the Gaussian planar diagrammatic expansion (Hermitian 
case in App. IB 1 bl non- Hermitian in Add. IB 2"b| . 



1. Hermitian random matrices 

a. Hermitian Green function 

Consider an N x N Hermitian random matrix H, with 
a corresponding probability measure d/x(H), in the large- 
TV limit [cf. (O] . One of the basic characteristics of this 
model is the "mean spectral density" (MSD), i.e., the dis- 
tribution of its (real) eigenvalues A^, averaged according 
to (...) = /(..• )dM(H), 



1 N 



(Bl) 



where S(x) is the real Dirac delta function. 

However, there exists an equivalent but handier object 
which comes in three versions: 

(i) The N x N "holomorphic Green function matrix," 



1 



Z H 



where Z = zl 



(B2) 



where z = x + iy is a complex argument. 

(ii) Its normalized trace is known as the "holomor- 
phic Green function" (other names: "resolvent," "Stjelt- 
jes transform" ) , 



N 



1 



A, 



(B3) 



(iii) In some situations an even more convenient quan- 
tity is the "holomorphic M-transform," 



M u (z) = zGu(z) - 1, 



(B4) 



whose expansion around z = oo (if it exists) is a power 
series of the (positive) "moments," 



Mh(z) = E^" 

n>l 
1 



m H ,n = ^Tr(H") = / dxp H (x)x r ' 

^ J cuts 



(B5a) 
(B5b) 



The relationship of the holomorphic Green function to 
the MSD is the following: 

(i) For any finite N, the former is a meromorphic func- 
tion, with poles at the mean eigenvalues; but when N 
grows to infinity, these poles merge into continuous in- 
tervals on the x-axis, distributed according to pn(x), 
and thus Gn(z) becomes holomorphic on the whole com- 
plex plane except these real "cuts." Then, the MSD 
can be used to calculate the holomorphic Green function 
through 



dxpu(x)- 



1 



(B6) 



(ii) On the other hand, since the real Dirac delta 
is known to have a representation (the "Sokhotsky- 
Weierstrass formula" ) , 



S(x) 



-— lim ( — - 

27T1 c->0+ \ X + 



1C 



(B7) 



one infers also an opposite link, 
pn{x) = -^r lim (G H (x + ie) - G u (x - ie)) , (B8) 

Z7T1 c->0+ 

i.e., the MSD is the jump of the holomorphic Green func- 
tion across the cuts. This is schematically illustrated in 
Fig. [2TJ 



complex plane 



real asris 



cuts of 
J, eigenvalues 



FIG. 21: For Hermitian random matrices, the real MSD is 
evaluated from the complex Green function of a complex ar- 
gument in the vicinity of the real axis 



b. Dyson-Schwinger equations and Gaussian unitary 
ensemble 

There exist various techniques of calculating the holo- 
morphic Green function of a given Hermitian random ma- 
trix model H. In the following, the method of "planar 
diagrammatic expansion" and "Dyson-Schwinger (DS) 
equations" [ixl HH will be introduced and applied as an 
example to solve a simple model, the "Gaussian unitary 
ensemble" (GUE), 



d/i(H) oc e 2c 



r TrH- 



dH, 



(B9) 



i.e., such that all Reify, ImTJ^ for i < j are in- 
dependent Gaussian of zero mean and variance a 2 /N 

moreover 



(for i 



j) or a 2 /(2N) (for i < j); 



dH = HidHuYl^dCReHijWlmHij) is the flat mea- 



To begin with, the holomorphic Green function matrix 
B2p should be expanded around z = oo, 



(BIO) 



G H (z) = Z- 1 + (Z^HZ^HZ -1 ) + 

+ (Z^HZ^HZ^HZ^HZ -1 ) + 
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where only the even moments have been kept due to the 
evenness of the GUE measure (|B9[) . 

Therefore, in order to sum up (|B10|) . one needs to cal- 
culate at each (even) order n the averages (n-point cor- 
relation functions) ([H] kl k . 2 [H] ks k 4 ■ ■ ■ [H] fc2 „ _ 1 k 2n ) • Ac- 
cording to the Wick theorem, this is done by making 
all the contractions, i.e., forming all the possible pair- 
ings of the H-factors, and replacing each of them by 
an appropriate 2-point correlation function (propagator), 
([H] fc fc [H]fc ofcg+1 ). The GUE propagator is inferred 
from (|B9j) to be 

2 

({HUn} kl ) = ^SuS jk . (Bii) 

This procedure is best explained using a pictorial lan- 
guage, 

• • = [ Z ~% = z 6 V> 

i 3 

= ([H]y[H]fei) = jj-SuSjk, 
i j kl 

since then at each order in (jBlOj) . written in terms of 
matrix indices, one obtains a set of alternating horizon- 
tal segments (Z -1 ) and double dots (H), where the pairs 
of double dots should be connected with double arcs 
(propagators) — this series for n = 0, 2, 4 reads 




i kjt2 j 



([ H ]fcife 2 [H]fe 3 fc 4 ) 



i kik-2 k$4 fcgfcg kyks j 

i 1 i 1 

( [H] k i k 2 [H] k 3 k 4 [H] fe5 ke [H] fe7 fes ) 




( [H] j fc 2 [H] fc3 fc4 [H] fe5 fc(3 [H] kr fc 8 ) 




(The fc-indices are summed over.) Hence, (i) at n = 0, 
there is just one Z _1 (one horizontal line); (ii) at n = 2, 
there are three Z _1 (three horizontal lines) and two H 
(two double dots), so there is one possible contraction 
into one propagator (one double arc); (iii) at n = 4, one 
finds five Z _1 and four H, so there are three possible 
contractions into two propagators (two double arcs); one 
of them yields a non-planar graph; (iv) etc. 

One might worry that all the above infinite series of 
graphs would have to be computed in order to evaluate 
(|B10|) . but in the large- N limit there exists a powerful 
alternative — infinite N is the 't Hooft limit in which only 
the planar diagrams contribute to the sum [92^ , and these 
can be effectively summed as follows: 

Consider the "one-line-irreducible" (1LI) diagrams 
(i.e., those that cannot by disconnected into two sub- 
graphs by cutting a single horizontal line), and let Sh(z) 
(the N x N "self-energy matrix") be their generating 
function, 




This subset is related to all the diagrams in a twofold 
way: 

(i) On one hand, any planar graph has necessarily a 
form of a number of 1LI graphs connected with horizontal 
segments, 
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+ ..., 



i.e., symbolically, 



GhW 



= Z- 1 +Z- 1 E H (2)Z- 1 + 
+ Z- 1 S H (z)Z- 1 S H (2)Z- 1 
1 

~ Z-Sh(z)' 



(B12) 



which is the "first Dyson-Schwinger equation." 

(ii) On the other hand, for Gaussian random matrices 
and in the planar limit, it is true that any 1LI diagram 
originates from some other diagram by adding to it an 
external double arc, 



i.e. 



[Sn(z)]i 



AT 

E 

1,^2 — 



[Gh(*)] 



[H] ifcl [H; 



(B13) 



which is the "second Dyson-Schwinger equation." For 
instance, inserting here the GUE propagator (IB11I) . 



Sgue(z) = a 2 G G VE(z)l 



AT- 



(B14) 



Analogously to the holomorphic Green function (|B3|) . one 
defines the "self-energy" as the normalized trace, 



£h(*) = ^TrS H (z), 



(B15) 



and for the GUE it reads Sgue(z) = ^Ggue(z)- 

The first (|B12p and second (|B13|) DS equations are 
a set of two matrix equations for G-r{z) and Sh(z) — 
this is the fundamental tool exploited in this paper. 
For example, for the GUE, plugging (|B14|) into (|B12|) . 
and taking the normalized trace of both sides gives a 
quadratic equation for the holomorphic Green function 
whose solution with the proper asymptotics at infinity 



[i.e., Gu(z — > oo) ~ 1/z, which is equivalent to the nor- 
malization of the MSD, f , dxpn(x) = 1, cf. (|B6])] is 



1 



GgUe(z) = 'Z — n ( Z - 



2<ry/z + 2a) , (B16) 



• where the square roots are principal. Applying (|B8j) . one 
arrives at the famous "Wigner semicircle law," 



Pgue(x) 




for x G [-2a, 2a] 
otherwise. 



(B17) 

This completes an exposition of the DS method in the 
Hermitian world; below it will be extended to the more 
involved non-Hermitian case. 



2. Non-Hermitian random matrices 

a. Non-Hermitian Green functions 

A fundamental difference between Hermitian and non- 
Hermitian random matrices which makes the latter (de- 
note an example by X) both richer and more difficult is 
that the eigenvalues Xi are generically complex instead of 
real, thus in the large- A limit coalesce not into real cuts 
but some two-dimensional "domains" T> on the z-planc 
with the MSD, 



1 N 

i=l 



6 {2) 



(B18) 



Nonholomorphic Green function. The above defini- 
tion is analogous to (|B1[) except that the Dirac delta 
function is now complex instead of real. Consequently, 
the Sokhotsky-Weierstrass formula (|B7[) can no longer be 
used as its representation, and therefore, the concept of 
the Green function (|B3j) needs to be reestablished accord- 
ing to another valid representation, 



z, z 



— lim 



— — — lim 



7T e-X) (U|2 _j_ £ 2\2 n e-K> \z\ 2 + 6 2 ' 

(B19) 

Inspired by this form, one introduces the "nonholomor- 
phic Green function" [93|-|97j], 



z - A,- 



Gx(z, z) = lim lim — / 

e _,o jv-kx, N ^ \ \ z _ \tf + e 2 I 

= lim lim — Tr/^ „. . N — — r - 

e^ON-yooN \ (zl N - X) (zl N - Xt) + £21^ 

(B20) 

(the matrix inversion will always be understood as 
A/B = AB _1 ), since then the MSD is simply propor- 
tional to its Wirtinger derivative with respect to z, 

p*(z,z) = ~-^Gx(z,z), for zeV. (B21) 

7T OZ 



N 
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Remark that, strictly speaking, the order of limits in 
(|B20I) is opposite to the right one: One should take e — > 
for any finite TV first, which can be directly performed 
everywhere except a finite number of points A^, but this 
would reduce the nonholomorphic Green function to the 
holomorphic one (|B3[) . Only for TV already infinite, and 
inside the (then continuous) mean spectral domain 2?, 
the limit e — > produces a nontrivial quantity. This 
intricacy will be disregarded in the following and an as- 
sumption made that for the matrix models investigated 
in this paper the limits do commute. 

Moreover, outside T>, the regulator e is unnecessary, 
and one is left with the holomorphic Green function, 



G x (z,z) = G x {z) 



for 



2?, 



(B22) 



which however — contrary to the Hermitian case — does 
not carry any information about the MSD. 

Alternatively to (|B20[) . it is sometimes more natural 
to use the "nonholomorphic A/-transform," defined anal- 
ogously to (|B4|) . 



Mx.{z,z) = zG x (z,z) - 1. 



(B23) 



Matrix-valued Green function. The nonholomorphic 
Green function (|B20|) encodes the MSD (jB2Tj) but still 
remains inconvenient to evaluate due to the denomina- 
tor quadratic in X unlike the linear denominator in the 
holomorphic version (IB3J . This problem ha s be en re- 
solved by "duplication" [98|, [§g| (cf. also [H [l0(3 for a 
slightly different approach), i.e., introducing a 27V x 27V 
matrix, 



Gy p1 *(z, z) = lim lim / — : — = 

X V ' e^O JV-Kx> \ 2 du P'- _ X du P'- 



(B24) 



where for short, 



^dupl. 



zl N 


ieljv 


ieljv 


zl N 



■j^dupl. 



X 




Oat 


xt 



(B25) 



This form exactly parallels the holomorphic one ()B2|) . 
with a denominator linear in the random matrix, for 
which the cost is the doubling of the matrix dimensions. 
The four TV x TV blocks of this and alike matrices will be 
distinguished by the following superscripts, 



r-idupl. 

U x = 



Gf 


Gf 


Gf 


Gf 



Taking the normalized "block-trace," 



bTr 




TrA 


TrB ' 


TrC 


TrD 



(B26) 



(B27) 



one obtains from (IB24[) the (2 x 2) "matrix-valued Green 
function," 



It is a fundamental object in the non-Hcrmitian world 
for two reasons: 

(i) Its upper left corner is precisely the desired non- 
holomorphic Green function (|B20j) . 



g^(z,z)=G x (z,z). 



(B29) 



In other words, the MSD has been encoded as (a deriva- 
tive of) a block of an object analogous to the holomorphic 
Green function (|B3|) . thus opening the door for meth- 
ods from Hermitian RMT (such as planar diagrams and 
DS equations; cf. Sec. IB 2 bp to be applied in the non- 
Hermitian case as well. [The lower right block is simply 
the complex conjugate of (|B29[) ; it does not bring any 
new information.] 

(ii) Furthermore, the off-diagonal blocks of (|B28[) . 



1 



lim lim — Tr 



e^ON^ooTV \{zl N -Xt) (zl N -X) + e 2 l N / ' 

(B30a) 

Gl z (z,z) = 

lim lim — Tr 



e ^ow^ooTV \(zl N -X)(zl N -X^) + e 2 l N / ' 

(B30b) 

which are equal and purely imaginary, so it is convenient 
to use their negated product, 

B x (z,z) = -g£(z,z)Gl z (z,z) G R+ U {0}, (B31) 

which has a twofold meaning: Firstly, it speaks about 
the cor relation between the left and right eigenvectors of 
X |l0lj . a property not investigated in this article. Sec- 
ondly, it is an "order parameter" : It vanishes outside the 
mean spectral domain T>, since the regulator is unneces- 
sary there, e = 0, and strictly positive inside this domain. 
Hence, once (|B31[) has been calculated inside T> (i.e., in 
the nonholomorphic sector), then 



B*(z,z) = 



(B32) 



is the equation of the borderline of T> in the Cartesian 
coordinates (x, y), i.e., z = x + iy € dT>. 

Quaternion Green function. Although not used in this 
paper, remark for completeness that the matrix- valued 
Green function [ (|B24p ~ (|B28|) ] may be naturally general- 
ized by replacing the infinitesimal parameter e with a 
complex (actually, real nonnegative is sufficient) argu- 
ment b (plus z is renamed a), 



gx(g) = ^^ Q&ljy 1 _ xdupL 
where the argument is a quaternion, 



(B33) 







<3H 





-bTrG^'P'- 

TV x 



(B28) 




a, b S 



(B34) 
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and the tensor product means that each of its four blocks 
is multiplied by 1/y. It is known as the "quaternion 
Green function" [102| . Il03j (cf. [H, [36| for its recent ap- 
plication). This extension allows to define, in analogy to 
(|A6|) . the functional inverse (in the quaternion space) of 
(IB33I) (the "quaternion Blue function"), which in turn 
leads to an extension of the free probability addition law 
(|A7j) from the Hcrmitian to non-Hcrmitian world. It also 
shows that just as in the Hermitian real MSD re- 

quires a complex Green function in the vicinity of the real 
axis (cf. Fig. l2"Tj) . so in the non- Hermitian case, a quater- 
nion Green function in the vicinity of the complex plane 
(cf. Fig. |2"2"|) needs to be employed in order to reproduce 
a complex density. 



quaternion space 




complex plane 



FIG. 22: For non-Hermitian random matrices, the complex 
MSD is evaluated from the quaternion Green function of a 
quaternion argument in the vicinity of the complex plane, i.e., 
from the matrix-valued Green function [ (|B24|1 - (|B26[1 . (|B28)> ]. 



where a bar over an index or its lack determines the block 
to which the index belongs, according to the convention 
(|B26|) (so here the nonzero propagator is between the 
upper left and lower right blocks). 

This is the necessary ingredient of the DS equations. 
The self-energy matrix is also duplicated, 



^dupl. — 



The hrst DS equation is analogous to ()B12|) . 

1 



Qjdu.pl . _ 



2Jdupl. _ ^Tjdupl. 



(B38) 



(B39) 



where for short, Z du P L = Z^ L (jB"25)) (one may set e = 
here as it will be seen that the DS equations by them- 
selves regularize the singularities inside the mean spectral 
domain). The second DS equation is analogous to (|B13|) 



N 



["y^dupl." 



i.e., 



^yidupl." 
^yidupl." 

^dupl. _ 



fel ,/C2 = l 



koj 



= [s dup1 -] 


-. = 0, 


Ojv 


<j 2 g zz i N 


a 2 g zz i N 


On 



(B40a) 
(B40b) 
(B40c) 



(B4I) 



Dyson-Schwinger equations and Ginibre unitary 
ensemble 



In order to get better acquainted with the DS tech- 
nique in the non-Hermitian case, the steps from Sec. IB 1 bl 
will now be applied to an exa mple, the (square) "Ginibre 
unitary ensemble" (GinUE) (T04m07j , 



d^X) oc e-5 T '( xtx )dX, 



(B35) 



i.e., all the ReXy, Im^ arc IID Gaussian random num- 
bers of zero mean and variance a 2 /(2N), for some a. 

One starts from writing down the nonzero propagators, 
immediately stemming from (|B35|) . 



<pq«[xt] w ) 



N 



o~u8jk- 



(B36) 



However, here one needs the propagator of the duplicated 
random matrix (|B25I) . 



<[X 



dupl." 



,[x 



dupl. 



(B37) 



Inserting ()B41|) into (|B39p . one finds a 2 x 2 matrix 
equation for the elements of the matrix-valued Green 
function of the GinUE, 



z 


-o 2 g 


ZZ 


-1 




-a 2 g zz 


z 


) 


■ 




i 








a 2 g zz 


2 _ a igzzgzz 


I a 2 g zz 


z 



(B42) 

Firstly, the equation for the order parameter 
B = Bx(z,z) (|B3ip reads therefore, 



B = 



o±B 



\z\ 2 + a A B) 



2 ' 



(B43) 



which has the "holomorphic solution" B = 0, valid out- 
side of D, and a relevant "nonholomorphic solution," 



B = \(l- V 



(B44) 
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It implies (|B32|) that the mean spectral domain T> for 
GinUE is the centered circle of radius 



\z\=a. 

Secondly, the upper left block of (|B42|) yields 



4r, for \z\ < cr, 
i for \z\ > a. 



\z\ 2 + a 4 B 



(B45) 



(B46) 



Hence, one recognizes the (trivial) holomorphic Green 
function l/z outside T>, and the nonholomorphic one, 
Gx(z,z) =~z/<j 2 , inside T>, which leads (|B21|) to the con- 
stant MSD, 



px.{z,z) 



-iff, for \z\ < cr, 
0, for \z\ > cr. 



(B47) 



This very method is used in this paper to solve much 
more complicated models. 



Appendix C: Mean spectral density of the 
equal-time covariance estimator for Toy Models 1, 
2a, 3 

This article is devoted mainly to the TLCE, but an 
important conclusion may be drawn from comparing the 
MSD of the ETCE and TLCE for a given model which 
is that the former is a better tool for discerning true 
spatial covariances, while the latter for temporal corre- 
lations. Therefore, even though the MSD of the ETCE 
is already known for all the models considered here, it is 
helpful to present these results and their derivation for 
three simplest examples, i.e., Toy Models 1, 2a and 3, for 
Gaussian assets. The pertinent master equation in all 
these cases is (l48l). 



MC histograms (solid)/theory (dashed), 
ttiterations=1000, ttbins=100, 
K=5, cr\=\, <x;=2, 0-5=4, 0-5=6, erf =8, p 1,2,3,4.5 =0.2 




r=0.02 (N=20, T=1000) 
r=0.1 (N=100, T=1000) 
r=0.5 (N=500, T=1000) 



20 



FIG. 23: Monte Carlo eigenvalues versus theory for the ETCE 
for Toy Model 2a. The graphs in the main figures are hardly 
distinguishable, hence the insets show them separately. 



which is the celebrated Marcenko-Pastur distribution 
(cf. Sec. II A 21 for its early appearance in econo- 
physics). [The first term in (|C3[) represents the zero 
modes, which appear only for r > 1, and which 
are derived by expanding (|C2[) around z — 0, which 
yields G = (1 - l/r)Q(r - l)(l/z) + O(z ), to which 



(|B7p needs to be applied. The normalization of the MSD 
essentially requires this contribution.] 



2. Toy Model 2a 



For Toy Model 2a (fT33|) (cf. Sec. IIIIC1|) , the master 
equation furnished with () 13411 and JVa(z) = 1 + l/z be- 



1. Toy Model 1 



For Toy Model 1 (fTT8|) (cf. Scc. ImB) . the master equa- 
tion becomes quadratic for G = G c (z), 

ra 2 zG 2 + (-z + a 2 (l- r)) G + 1 = 0, (CI) 

whose solution with the proper asymptotics at infinity 
[cf. the discussion before eq. (|B16|) ] reads 



G 



1 



2rer 2 z 



(2 - cr 2 (l - r )-y/z- x+^z - x-) , (C2) 



where x± = cr 2 (l ± \/t) 2 , and where the square roots are 
principal. Hence, the MSD follows from (|B8|) to be 



Pc(x) 



(1 - l/r)9(r- l)5(x)+ 



\f (X-L— X)(X— X-) r r -I 

— -\ forxG [x-,x+ , 

2-nr<7 z x ' L ■? -t- J ^ 

0, otherwise, 



(C3) 



M 



K 

E 

k=l rM+l 



(C4) 



which is a polynomial equation of order (K + 1) for the 
holomorphic M-transform M = M c (z), whose solution — 
with [pi]). dHJ] applied to it— yields the desired MSD. 

A more detailed analysis of this model has been done 
elsewhere (cf. e.g. [27|). but in Fig. [23] histograms of 
Monte Carlo eigenvalues are compared to the solutions 
of ([C4|, for K = 5 with a\ = 1, a\ = 2, a\ = 4, a\ = 6, 
erf = 8 and equal relative multiplicities 0.2, as well as 
three values of r = 0.02,0.1,0.5, discovering perfect 
agreement between the two. One recognizes how de- 
creasing r makes the peaks around the true variances 
narrower, allowing their determination from this MSD. 
This should be contrasted with Fig. [6] (c), which shows 
that even for very small r the separate variances are not 
visible in the spectrum of the TLCE. 
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FIG. 24: Monte Carlo eigenvalues versus theory for the ETCE 
for Toy Model 3. 



3. Toy Model 3 



a. Mean spectral density 



For Toy Model 3 (fl"42j) (cf. Scc lmD]) . the master equa- 
tion with Nc(z) = 1 + 1/ z and 



M A (z) 



1/2 



dp- 



1 



- 1 



(C5) 



± 



tanh 77- a ; 

2r V o- z 



COth 77- 



b. Edges of the spectrum (new result) 



It is a general principle [1081 . 1 1091 ] that for any Hermi- 
tian random matrix model H, the edges a;* of its mean 
spectrum (in the thermodynamic limit) are given by the 
divergences of the Green function, dGn(z) /dz\ z — Xt = oo. 
Combining this formula with (|C6j) . one arrives at a sextic 
polynomial equation for the positions of the edges, 



( X 2 - 1) xl, - 6r X {x 2 - 1) 4*+ 
+ (3(l-r 2 )-(2 + 9r 2 ) X 2 + 12rV)4.+ 
+ 2r X (i (1 + 2r 2 ) - (5 + 2r 2 ) X 2 - 4r 2 X 4 )4*+ 
+ ( - 3 (l + 7r 2 + r 4 ) + (l + 26r 2 - 9r 4 ) X 2 + 
+ r 2 (l + 12r 2 ) x 4 )^L" 

- 2r X (3 (1 - r 2 ) (2 + r 2 ) + r 2 (5 + 3r 2 ) x 2 )^*+ 
+ (l-r 2 ) 2 (l-r 2 +r 2 X 2 )=0, 



(C7) 

where real and nonnegative solutions should be picked, 
and it is a numerical observation (to be proven) that there 
are two of them, as expected. Figure [9] shows these edges 
as functions of r, compared with the analogous quantities 
for the TLCE. 



Appendix D: Numerical evaluation of the mean 
spectral density of the time-lagged estimator for Toy 
Model 3 



(the sign is irrelevant, the square roots are principal), is 
a quartic polynomial equation for G = G c (z), 

r 2 z 2 G 4 a ~2r Z<7 { X z a + r) G%+ 
+ (zl + 4rxz a + r 2 -l)Gl- (C6) 
- 2 (z a + r X ) G a + 1 = 0, 

where for short, x = coth(l/r), where the argument and 
the holomorphic Green function are rescaled by the vari- 
ance, z„ = zja 2 and G a = a 2 G. 

This result has been found before (28|, nevertheless 
Fig. [2H presents some tests of (|C6[) against Monte Carlo 
data, for r = 2.5 and three values of r = 0.02,0.1,0.5, 
with perfect concord. A crucial observation is that the 
MSD is "MP-shaped," i.e., it may be impossible to distin- 
guish from the case with no temporal correlations, which 
means that the ETCE is not an effective tool in discern- 
ing whether a given empirical set of data exhibits true 
temporal covariances or not. On the other hand, com- 
paring Fig. [1] with Figs. QT1[]21[I2] reveals that the spec- 
trum of the TLCE is highly sensitive to whether there 
are intrinsic temporal correlations or not, and of what 
strength. 



This appendix briefly describes a simple nu- 
merical algorithm used to calculate the MSD, 
Pc(t)( x, y) = k^G = ±(d x X - d y Y + i(d y X + d x Y)) 
(|B21|) of the TLCE for Toy Mode l 3 (fT42|) ( cf. Sec . [HTDj) 
from the master equations [ (|145a| . (|145bj) . p46ap - (fl47l) ]. 
The algorithm relies on certain observations, theorems 
and conjectures which are made for some special values 
of the parameters of the model (i.e., r — 0.5, a = 1, 
r = 5, t = 10), but are supposed to hold in general. 
Even though Toy Model 3 is very simple, the algorithm 
captures certain generic properties any numerical ap- 
proach to the MSD of more realistic financial models 
will probably share. 

The basic observation is that the argument z appears 
in the master equations only on the RHS of (|145b[) . hence 
the method will be to (i) pick a value of G (from a certain 
regular lattice), and solve (|145ap for h > 0, (ii) insert 
these G and h into (|145bj) . thereby obtaining the value 
of z corresponding to the assumed value of G. 

In this program, one is therefore faced with two prob- 
lems: (i) How to solve (|145a)l for h, with given Gl 
(Scc.[DTJ) (ii) How to numerically compute the deriva- 
tive dzG from the set of values z obtained from a regular 
lattice of the values of Gl (Sec. [D2j) 
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1. Solving the first master equation (|145ap for h 

An initial observation is that solving (|145a| for h needs 
to be done in some non-iterative way, because it seems 



to behave chaotically upon iterating. 

Before attempting a solution, one has to know how 
to calculate the integrals Fi t2 (G, ~h) [ (|146a| . (|146b|) ] for 
given arguments. For this, note that for h > the roots 
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FIG. 26: Toy Model 3: Interpolation of the MSD fSec.lD2)l. The refining divisions sgj. = s£S. = 5. Figs, (c) and (d) are zooms 
of (a) and (b), respectively, into [-0.25,0.25] x [-0.25,0.25]. 



of W(u) (|147p can never lie on the unit circle G(0, 1); it 
is clear that if uo was a root satisfying uq = 1/uo, then 
W(u ) = ul +1 (|1 - A 2 u + u 2 + rA l Gu t + l \ 2 + r 2 A\h) 
> 0, which is a contradiction. Hence, there is certainly 
no divergence in the contour integrals present in the 
master equations. Moreover, it also follows from (|147p 
that the roots of W(u) come in pairs, (u n , for 
n = 1, ...,£ + 1. These two properties mean that al- 
ways half of the roots of W(u) lie inside C(0, 1), and the 
other half outside; only the internal ones contribute to 
the contour integrals. Therefore, by expressing the in- 
tegrals through the residues at these internal roots, one 
acquires precise control over the values of Fi^{G, h), and 
obtains a set of algebraic (i.e., no longer integral) master 
equations, much more suitable for numerical evaluation. 



points G for which (|145a[) has at least one solution h > 
form a bounded area T>q. This is crucial from the point of 
view of numerical analysis, as one needs to search through 
only a restricted domain in the G-plane. (ii) For any 
G € T>g, there can only be either one or two solutions h. 
This is depicted in Fig. [25] (a), where the points are the 
subset of the initial G-lattice forming 2? G , and the two 
shades of gray correspond to one or two solutions. These 
solutions form two continuous branches; Figs. [25] [(b), (c)] 
show the G-domains of these two branches, while their 
plots are presented in Figs. [5S] (d) (both branches glued 
together to form a continuous surface %) and [(e), (f)] 
(the branches separately). 



Solving (|145aj) for h is thus performed as follows: Form 
a regular lattice in the G-plane (for example, a rectan- 



gle [X r 



J" 



divided into s 



Rc 



..Iiii 



stripes in each direction^ respectively) , and for each point 
G, scan an interval [0, /i max ] searching for solutions h of 
(|145ap . The upper bound of solutions /i mM , as well as 
the lattice parameters, are tuned in by trial and error. 

Performed this program numerically for the selected 
parameters of the model, two important observations are 
made, which are conjectured to hold in general: (i) The 



Finally, one may substitute each G from the lattice and 
its corresponding h (from each of the two branches) into 
(|145b[) . obtaining the value of z. The points z computed 
from both branches (reproducing theoretically the mean 
spectral domain T>) arc marked by proper shades of gray 
in Figs. [26] [(a), (c)], and shown on the background of the 
Monte Carlo eigenvalues of the TLCE. This describes a 
bijection H — > T>. 
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2. Interpolation and calculating the mean spectral 
density 

Once the above numerical analysis is completed, one 
still faces another obstacle: Given the facts that one typ- 
ically picks the initial G-lattice somewhat larger than the 
domain Vq plus that it is time-consuming to solve (|145al) 
at each lattice point — one generically obtains at the end 
too few points z to calculate from them the MSD of the 
TLCE with a satisfactory accuracy. 

For this reason, linear interpolation have been resorted 
to in order to create worthwhile MSD graphs: (i) For 
each of the two branches separately, select those lattice 



points from their G-domains which form the four cor- 
ners of an elementary rectangle of the lattice which en- 
tirely lies within the respective domain, (ii) Create a 
sub- lattice by dividing the rectangle in question into sf° t , 
s i"t ( = 5 here) stripes in each direction, (hi) Compute 
the value of z at each point of this sub-lattice by linear 
interpolation from the values of z at the four corners of 
the sub-lattice; this is presented in Figs. [26] [(b), (d)]. 
In addition to z, it is straightforward to calculate the 
derivatives dxx, dxy, dyx, dyy, hence also the inverse 
derivatives d x X, d y X, d x Y, d y Y through an appropri- 
ate Jacobian — which is equivalent to knowing the MSD 
Pc(t)(x,y), plotted in Fig.[H 
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